Analytic first and second derivatives of the energy in the fragment molecular orbital method combined with molecular mechanics

Analytic first and second derivatives of the energy are developed for the fragment molecular orbital method interfaced with molecular mechanics in the electrostatic embedding scheme at the level of Hartree-Fock and density functional theory. The importance of the orbital response terms is demonstrated. The role of electrostatic embedding upon molecular vibrations is analyzed, comparing force field and quantum mechanical treatments for an ionic liquid and a solvated protein. The method is applied for 100 protein conformations sampled in molecular dynamics (MD) to take into account the complexity of a flexible protein structure in solution, and a good agreement with experimental data is obtained: Frequencies from an experimental infrared (IR) spectrum are reproduced within 17 cm⁻¹ .     

Theoretical study of the prion protein based on the fragment molecular orbital method

We performed fragment molecular orbital (FMO) calculations to examine the molecular interactions between the prion protein (PrP) and GN8, which is a potential curative agent for prion diseases. This study has the following novel aspects: we introduced the counterpoise method into the FMO scheme to eliminate the basis set superposition error and examined the influence of geometrical fluctuation on the interaction energies, thereby enabling rigorous analysis of the molecular interaction between PrP and GN8. This analysis could provide information on key amino acid residues of PrP as well as key units of GN8 involved in the molecular interaction between the two molecules. The present FMO calculations were performed using an original program developed in our laboratory, called “Parallelized ab initio calculation system based on FMO (PAICS)”. © 2009 Wiley Periodicals, Inc. J Comput Chem 2009     

CH/pi hydrogen bonds determine the selectivity of the Src homology 2 domain to tyrosine phosphotyrosyl peptides: an ab initio fragment molecular orbital study

The CH/pi hydrogen bond is a weak molecular force occurring between CH groups (soft acids) and pi-systems (soft bases), and has been recognized to be important in the interaction of proteins with their specific ligands. For instance, it is well known that Src homology-2 protein (SH2) recognizes its specific pTyr peptide in two key regions, pTyr-binding region and specificity-determining region, by the use of attractive molecular forces, including the CH/pi hydrogen bond. We hypothesized that the CH/pi hydrogen bond plays a key role in determining the selectivity of SH2 proteins, and studied this issue by the ab initio fragment molecular orbital (FMO) method. The FMO calculations were carried out, at the HF/6-31G* and MP2/6-31G* level, for SH2 domains of Src, Grb2, P85alpha(N), Syk, and SAP, in complex with corresponding pTyr peptides. CH/pi hydrogen bonds have in fact been found to be important in stabilizing the structure of the complexes. We conclude that the CH/pi hydrogen bond plays an indispensable role in the recognition of SH2 domains with their specific pTyr peptides, thus playing a vital role in the signal transduction system.     

Ab initio quantum chemical calculation of electron density,electrostatic potential,and electric field of biomolecule based on fragment molecular orbital method

Efficient quantum chemical calculations of electrostatic properties, namely, the electron density (EDN), electrostatic potential (ESP), and electric field (EFL), were performed using the fragment molecular orbital (FMO) method. The numerical errors associated with the FMO scheme were examined at the HF, MP2, and RI‐MP2 levels of theory using 4 small peptides. As a result, the FMO errors in the EDN, ESP, and EFL were significantly smaller than the magnitude of the electron correlation effects, which indicated that the FMO method provides sufficiently accurate values of electrostatic properties. In addition, an attempt to reduce the computational effort was proposed by combining the FMO scheme and a point charge approximation. The error due to this approximation was examined using 2 proteins, prion protein and human immunodeficiency virus type 1 protease. As illustrative examples, the ESP values at the molecular surface of these proteins were calculated at the MP2 level of theory.     

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.     

A new hierarchical parallelization scheme: generalized distributed data interface (GDDI), and an application to the fragment molecular orbital method (FMO)

A two-level hierarchical scheme, generalized distributed data interface (GDDI), implemented into GAMESS is presented. Parallelization is accomplished first at the upper level by assigning computational tasks to groups. Then each group does parallelization at the lower level, by dividing its task into smaller work loads. The types of computations that can be used with this scheme are limited to those for which nearly independent tasks and subtasks can be assigned. Typical examples implemented, tested, and analyzed in this work are numeric derivatives and the fragment molecular orbital method (FMO) that is used to compute large molecules quantum mechanically by dividing them into fragments. Numeric derivatives can be used for algorithms based on them, such as geometry optimizations, saddle-point searches, frequency analyses, etc. This new hierarchical scheme is found to be a flexible tool easily utilizing network topology and delivering excellent performance even on slow networks. In one of the typical tests, on 16 nodes the scalability of GDDI is 1.7 times better than that of the standard parallelization scheme DDI and on 128 nodes GDDI is 93 times faster than DDI (on a multihub Fast Ethernet network). FMO delivered scalability of 80-90% on 128 nodes, depending on the molecular system (water clusters and a protein). A numerical gradient calculation for a water cluster achieved a scalability of 70% on 128 nodes. It is expected that GDDI will become a preferred tool on massively parallel computers for appropriate computational tasks.     

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     

Maximal orbital analysis of molecular wavefunctions

We describe a new way to decompose one-electron orbitals of a molecule into atom-centered or fragment-centered orbitals by an approach that we call “maximal orbital analysis” (MOA). The MOA analysis is based on the corresponding orbital transformation (COT) that has the unique mathematical property of maximizing any sub-trace of the overlap matrix, in Hilbert metric sense, between two sets of nonorthogonal orbitals. Here, one set comprises the molecule orbitals (Hartree–Fock, Kohn–Sham, complete-active-space, or any set of orthonormal molecular orbitals), the other set comprises the basis functions associated with an atom or a group of atoms. We show in prototypical molecular systems such as a water dimer, metal carbonyl complexes, and a mixed-valent transition metal complex, that the MOA orbitals capture very well key aspects of wavefunctions and the ensuing chemical concepts that govern electronic interactions in molecules. © 2019 Wiley Periodicals, Inc.     

Studies in the paired orbital method

The energy expression and its derivatives have been derived in the paired orbital (PO) method. Several relations have been derived for the values of the basic functions at the Hartree-Fock limit. These relations permit to draw conclusions about the behavior of the energy expression as a function of the nonlinear parameters around this limit point and it shows that one can expect an improvement in energy as compared to the Hartree-Fock value using the nonlinear parameters as variational parameters.Dedicated to Professor J. Koutecký on the occasion of his 65th birthday     

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.     

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.     

Generalized hybrid orbital for the treatment of boundary atoms in combined quantum mechanical and molecular mechanical calculations using the semiempirical parameterized model 3 method

The application of combined quantum mechanical (QM) and molecular mechanical methods to large molecular systems requires an adequate treatment of the boundary between the two approaches. In this article, we extend the generalized hybrid orbital (GHO) method to the semiempirical parameterized model 3 (PM3) Hamiltonian combined with the CHARMM force field. The GHO method makes use of four hybrid orbitals, one of which is included in the QM region in self-consistent field optimization and three are treated as auxiliary orbitals that do not participate in the QM optimization, but they provide an effective electric field for interactions. An important feature of the GHO method is that the semiempirical parameters for the boundary atom are transferable, and these parameters have been developed for a carbon boundary atom consistent with the PM3 model. The combined GHO-PM3/CHARMM model has been tested on molecular geometry and proton affinity for a series of organic compounds.Acknowledgement We thank the National Institutes of Health for support of this research.Contribution to the Jacopo Tomasi Honorary Issue     

Energy gradients in combined fragment molecular orbital and polarizable continuum model (FMO/PCM) calculation

The analytic energy gradients for the combined fragment molecular orbital and polarizable continuum model (FMO/PCM) method are derived and implemented. Applications of FMO/PCM geometry optimization to polyalanine show that the structures obtained with the FMO/PCM method are very close to those obtained with the corresponding full ab initio PCM methods. FMO/PCM (RHF/6‐31G* level) is used to optimize the solution structure of the 304‐atom Trp‐cage miniprotein and the result is in agreement with NMR experiments. The key factors determining the relative stability of the α‐helix, β‐turn and the extended form in solution are elucidated for polyalanine. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2010     

Full variational molecular orbital method: Application to the positron-molecule complexes

Optimal Gaussian-type orbital (GTO) basis sets of positron and electron in positron-molecule complexes are proposed by using the full variational treatment of molecular orbital (FVMO) method. The analytical expression for the energy gradient with respect to parameters of positronic and electronic GTO such as the orbital exponents, the orbital centers, and the linear combination of atomic orbital (LCAO) coefficients, is derived. Wave functions obtained by the FVMO method include the effect of electronic or positronic orbital relaxation explicitly and satisfy the virial and Hellmann–Feynman theorems completely. We have demonstrated the optimization of each orbital exponent in various positron-atomic and anion systems, and estimated the positron affinity (PA) as the difference between their energies. Our PA obtained with small basis set is in good agreement with the numerical Hartree–Fock result. We have calculated the OH⁻ and [OH⁻; e⁺] species as the positron-molecular system by the FVMO method. This result shows that the positronic basis set not only becomes more diffuse but also moves toward the oxygen atom. Moreover, we have applied this method to determine both the nuclear and electronic wave functions of LiH and LiD molecules simultaneously, and obtained the isotopic effect directly. © 1998 John Wiley & Sons, Inc. Int J Quant Chem 70: 491–501, 1998     

Theoretical analysis of binding specificity of influenza viral hemagglutinin to avian and human receptors based on the fragment molecular orbital method

The hemagglutinin (HA) protein of the influenza virus binds to the host cell receptor in the early stage of viral infection. A change in binding specificity from avian 2-3 to human 2-6 receptor is essential for optimal human-to-human transmission and pandemics. Therefore, it is important to reveal the key factors governing the binding affinity of HA-receptor complex at the molecular level for the understanding and prediction of influenza pandemics. In this work, on the basis of ab initio fragment molecular orbital (FMO) method, we have carried out the interaction energy analysis of HA-receptor complexes to quantitatively elucidate the binding specificity of HAs to avian and human receptors. To discuss the binding property of influenza HA comprehensively, a number of HAs from human H1, swine H1, avian H3 and avian H5 viruses were analyzed. We performed detailed investigations about the interaction patterns of complexes of various HAs and receptor analogues, and revealed that intra-molecular interactions between conserved residues in HA play an important role for HA-receptor binding. These results may provide a hint to understand the role of conserved acidic residues at the receptor binding site which are destabilized by the electrostatic repulsion with sialic acid. The calculated binding energies and interaction patterns between receptor and HAs are consistent with the binding specificities of each HA and thus explain the receptor binding mechanism. The calculated results in the present analysis have provided a number of viewpoints regarding the models for the HA-receptor binding specificity associated with mutated residues. Examples include the role of Glu190 and Gln226 for the binding specificity of H5 HA. Since H5 HA has not yet been adapted to human receptor and the mechanism of the specificity change is unknown, this result is helpful for the prediction of the change in receptor specificity associated with forthcoming possible pandemics.     

Fragment molecular orbital study of the electronic excitations in the photosynthetic reaction center of Blastochloris viridis

All electron calculations were performed on the photosynthetic reaction center of Blastochloris viridis, using the fragment molecular orbital (FMO) method. The protein complex of 20,581 atoms and 77,754 electrons was divided into 1398 fragments, and the two‐body expansion of FMO/6‐31G* was applied to calculate the ground state. The excited electronic states of the embedded electron transfer system were separately calculated by the configuration interaction singles approach with the multilayer FMO method. Despite the structural symmetry of the system, asymmetric excitation energies were observed, especially on the bacteriopheophytin molecules. The asymmetry was attributed to electrostatic interaction with the surrounding proteins, in which the cytoplasmic side plays a major role. © 2009 Wiley Periodicals, Inc. J Comput Chem 2010     

Parallelized integral-direct CIS(D) calculations with multilayer fragment molecular orbital scheme

We have developed a parallelized integral-direct code of the perturbative doubles correction for configuration interaction with singles, proposed as CIS(D) by Head-Gordon et al. (Chem Phys Lett 219:21, 1994). The CIS(D) method provides the energy corrections both of the relaxation and differential correlation for the respective CIS excited states. The implementation of CIS(D) is based on our original algorithm for the second-order Møller–Plesset perturbation (MP2) calculations (Mochizuki et al. in Theor Chem Acc 112:442, 2004). There is no need to communicate bulky intermediate data among worker processes of the parallelized execution. This CIS(D) code is then incorporated into a developer version of ABINIT-MP program, in order to improve the overestimation in excitation energies calculated by the CIS method in conjunction with the multilayer fragment molecular orbital scheme (MLFMO-CIS) (Mochizuki et al. in Chem Phys Lett 406:283, 2005). The MLFMO-CIS(D) method is first used in evaluating the lowest n\(\pi^{*}\) excitation energy of the hydrated formaldehyde. The photoactive yellow protein (PYP) is the second target of MLFMO-CIS(D) calculation. Through these applications, it is shown that the CIS(D) correction improves the CIS results favorably.     

Ab initio quantum mechanical study of the binding energies of human estrogen receptor alpha with its ligands: an application of fragment molecular orbital method

We have theoretically examined the relative binding affinities (RBA) of typical ligands, 17beta-estradiol (EST), 17alpha-estradiol (ESTA), genistein (GEN), raloxifene (RAL), 4-hydroxytamoxifen (OHT), tamoxifen (TAM), clomifene (CLO), 4-hydroxyclomifene (OHC), diethylstilbestrol (DES), bisphenol A (BISA), and bisphenol F (BISF), to the alpha-subtype of the human estrogen receptor ligand-binding domain (hERalpha LBD), by calculating their binding energies. The ab initio fragment molecular orbital (FMO) method, which we have recently proposed for the calculations of macromolecules such as proteins, was applied at the HF/STO-3G level. The receptor protein was primarily modeled by 50 amino acid residues surrounding the ligand. The number of atoms in these model complexes is about 850, including hydrogen atoms. For the complexes with EST, RAL, OHT, and DES, the binding energies were calculated again with the entire ERalphaLBD consisting of 241 residues or about 4000 atoms. No significant difference was found in the calculated binding energies between the model and the real protein complexes. This indicates that the binding between the protein and its ligands is well characterized by the model protein with the 50 residues. The calculated binding energies relative to EST were very well correlated with the experimental RBA (the correlation coefficient r=0.837) for the ligands studied in this work. We also found that the charge transfer between ER and ligands is significant on ER-ligand binding. To our knowledge, this is the first achievement of ab initio quantum mechanical calculations of large molecules such as the entire ERalphaLBD protein.