Recent developments in QM/MM methods towards open-boundary multi-scale simulations

Combined quantum mechanics/molecular mechanics (QM/MM) methods have been widely used in multi-scale modelling and simulations of physical, chemical and biological processes in complex environments. In this review, we provide an overview of the recently developed QM/MM algorithms, with emphasis on our works, towards the ultimate goal of establishing an open boundary between the QM and MM subsystems. The open boundary is characterised by on-the-fly exchanges of partial charges and atoms between the QM and MM subsystems, allowing us to focus on the small QM subsystem of primary interest in dynamics simulations. An open-boundary scheme has the promise to the utilisations of small QM subsystems, high-levels of QM theory and long simulation times, which can potentially lead to new insights.     

Combined quantum and molecular mechanics calculations on metalloproteins

The combination of quantum mechanics and molecular mechanics (QM/MM) methods is one of the most promising approaches to study the structure, function and properties of proteins. The number of QM/MM applications on metalloproteins is steadily increasing, especially studies with density functional methods on redox-active metal centres. Recent developments include new parameterised methods to treat covalent bonds between the quantum and classical systems, methods to obtain free energy from QM/MM results, and the combination of quantum chemistry and protein crystallography.     

Multiscale methods for macromolecular simulations

In this article we review the key modeling tools available for simulating biomolecular systems. We consider recent developments and representative applications of mixed quantum mechanics/molecular mechanics (QM/MM), elastic network models (ENMs), coarse-grained molecular dynamics, and grid-based tools for calculating interactions between essentially rigid protein assemblies. We consider how the different length scales can be coupled, both in a sequential fashion (e.g. a coarse-grained or grid model using parameterization from MD simulations), and via concurrent approaches, where the calculations are performed together and together control the progression of the simulation. We suggest how the concurrent coupling approach familiar in the context of QM/MM calculations can be generalized, and describe how this has been done in the CHARMM macromolecular simulation package.     

The non-additive contribution of hydroxyl substituents to Akt kinase–apigenin affinity

The natural product apigenin is a flavonoid derivative substituted by three hydroxyl functional groups at positions 4′, 5 and 7 [OH(4′), OH(5) and OH(7)] of the basic flavonoid skeleton, which has shown strong inhibition on the development, proliferation and invasion of tumour cells by binding specifically to Akt kinase to inactivate the Akt signalling pathway. In this study, a typical non-additivity of the three hydroxyl substituents’ contributions to Akt–apigenin binding affinity is demonstrated by combination of four empirical scoring functions, molecular dynamics simulations, molecular mechanics-Poisson–Boltzmann/surface area (MM/PBSA) analyses and hybrid quantum mechanics/molecular mechanics (QM/MM) calculations. It is found that (i) the empirical scoring functions are incapable of properly reflecting the non-additivity feature, which, however, can be well described by the more rigorous MM/PBSA and QM/MM methods, (ii) the hydroxyl group contributions to ligand binding affinity are deviated significantly from linear additive model due to the strong conjugate effect and σ-effect among them, that is, the co-contribution of the three hydroxyl groups is far less than the sum of their individual contributions and (iii) as might be expected, a strong interactive effect is observed for the two adjacent substituents OH(5) and OH(7) as compared with that of distant OH(5) and OH(4′) as well as OH(7) and OH(4′). In addition, the structural basis, energetic property and molecular mechanism of the non-additivity feature are also explored in detail using the natural population analysis and quantum mechanical calculations.     

Binding free energy calculation with QM/MM hybrid methods for Abl-Kinase inhibitor

We report a Quantum mechanics/Molecular Mechanics–Poisson-Boltzmann/ Surface Area (QM/MM-PB/SA) method to calculate the binding free energy of c-Abl human tyrosine kinase by combining the QM and MM principles where the ligand is treated quantum mechanically and the rest of the receptor by classical molecular mechanics. To study the role of entropy and the flexibility of the protein ligand complex in a solvated environment, molecular dynamics calculations are performed using a hybrid QM/MM approach. This work shows that the results of the QM/MM approach are strongly correlated with the binding affinity. The QM/MM interaction energy in our reported study confirms the importance of electronic and polarization contributions, which are often neglected in classical MM-PB/SA calculations. Moreover, a comparison of semi-empirical methods like DFTB-SCC, PM3, MNDO, MNDO-PDDG, and PDDG-PM3 is also performed. The results of the study show that the implementation of a DFTB-SCC semi-empirical Hamiltonian that is derived from DFT gives better results than other methods. We have performed such studies using the AMBER molecular dynamic package for the first time. The calculated binding free energy is also in agreement with the experimentally determined binding affinity for c-Abl tyrosine kinase complex with Imatinib.     

QM/MM model study on properties and structure of some antibiotics in gas phase: Comparison of energy and NMR chemical shift

The combination of Quantum Mechanics (QM) and Molecular Mechanics (MM) methods has become an alternative tool for many applications for which pure QM and MM are not suitable. The QM/MM method has been used for different types of problems, for example: structural biology, surface phenomena, and liquid phase. In this paper, we have used these methods for antibiotics and then we compare results. The calculations were done by the full ab initio method (HF/3-21G) and the (HF/STO-3G) and QM/MM (ONIOM) method with HF (3-21G)/AM1/UFF and HF (STO-3G)/AM1/UFF. We found the geometry that has obtained by the QM/MM method to be very accurate, and we can use this rapid method in place of time consuming ab initio methods for large molecules. Comparison of energy values in the QM/MM and QM methods is given. In the present work, we compare chemical shifts and conclude that the QM/MM method is a perturbed full QM method. The work has been done on penicillin, streptomycin, benzyl penicillin, neomycin, kanamycin, gentamicin, and amoxicillin.     

Determinants of reactivity and selectivity in soluble epoxide hydrolase from quantum mechanics/molecular mechanics modeling

Soluble epoxide hydrolase (sEH) is an enzyme involved in drug metabolism that catalyzes the hydrolysis of epoxides to form their corresponding diols. sEH has a broad substrate range and shows high regio- and enantioselectivity for nucleophilic ring opening by Asp333. Epoxide hydrolases therefore have potential synthetic applications. We have used combined quantum mechanics/molecular mechanics (QM/MM) umbrella sampling molecular dynamics (MD) simulations (at the AM1/CHARMM22 level) and high-level ab initio (SCS-MP2) QM/MM calculations to analyze the reactions, and determinants of selectivity, for two substrates: trans-stilbene oxide (t-SO) and trans-diphenylpropene oxide (t-DPPO). The calculated free energy barriers from the QM/MM (AM1/CHARMM22) umbrella sampling MD simulations show a lower barrier for phenyl attack in t-DPPO, compared with that for benzylic attack, in agreement with experiment. Activation barriers in agreement with experimental rate constants are obtained only with the highest level of QM theory (SCS-MP2) used. Our results show that the selectivity of the ring-opening reaction is influenced by several factors, including proximity to the nucleophile, electronic stabilization of the transition state, and hydrogen bonding to two active site tyrosine residues. The protonation state of His523 during nucleophilic attack has also been investigated, and our results show that the protonated form is most consistent with experimental findings. The work presented here illustrates how determinants of selectivity can be identified from QM/MM simulations. These insights may also provide useful information for the design of novel catalysts for use in the synthesis of enantiopure compounds.     

QM/MM model study on properties and structure of some antibiotics in gas phase: Comparison of energy and NMR chemical shift

The combination of Quantum Mechanics (QM) and Molecular Mechanics (MM) methods has become an alternative tool for many applications for which pure QM and MM are not suitable. The QM/MM method has been used for different types of problems, for example: structural biology, surface phenomena, and liquid phase. In this paper, we have used these methods for antibiotics and then we compared results. The calculations were done by the full ab initio method (HF/3-21G) and the (HF/STO-3G) and QM/MM (ONIOM) method with HF (3-21G)/AM1/UFF and HF (STO-3G)/AM1/UFF. We found the geometry obtained by the QM/MM method to be very accurate, and we can use this rapid method in place of time consuming ab initio methods for large molecules. Comparison of energy values in the QM/MM and QM methods is given. In the present work, we compare chemical shifts and conclude that the QM/MM method is a perturbed full QM method. The work has been done on penicillin, streptomycin, benzyl penicillin, neomycin, kanamycin, gentamicin, and amoxicillin.     

Harder,better, faster,stronger: Large-scale QM and QM/MM for predictive modeling in enzymes and proteins

Computational prediction of enzyme mechanism and protein function requires accurate physics-based models and suitable sampling. We discuss recent advances in large-scale quantum mechanical (QM) modeling of biochemical systems that have reduced the cost of high-accuracy models. Tradeoffs between sampling and accuracy have motivated modeling with molecular mechanics (MM) in a multiscale QM/MM or iterative approach. Limitations to both conventional density-functional theory and classical MM force fields remain for describing noncovalent interactions in comparison to experiment or wavefunction theory. Because predictions of enzyme action (i.e. electrostatics), free energy barriers, and mechanisms are sensitive to the protocol and embedding method in QM/MM, convergence tests and systematic methods for quantifying QM-level interactions are a needed, active area of development.     

A comparative study of trypsin specificity based on QM/MM molecular dynamics simulation and QM/MM GBSA calculation

Hydrogen bonding and polar interactions play a key role in identification of protein-inhibitor binding specificity. Quantum mechanics/molecular mechanics molecular dynamics (QM/MM MD) simulations combined with DFT and semi-empirical Hamiltonian (AM1d, RM1, PM3, and PM6) methods were performed to study the hydrogen bonding and polar interactions of two inhibitors BEN and BEN1 with trypsin. The results show that the accuracy of treating the hydrogen bonding and polar interactions using QM/MM MD simulation of PM6 can reach the one obtained by the DFT QM/MM MD simulation. Quantum mechanics/molecular mechanics generalized Born surface area (QM/MM-GBSA) method was applied to calculate binding affinities of inhibitors to trypsin and the results suggest that the accuracy of binding affinity prediction can be significantly affected by the accurate treatment of the hydrogen bonding and polar interactions. In addition, the calculated results also reveal the binding specificity of trypsin: (1) the amidinium groups of two inhibitors generate favorable salt bridge interaction with Asp189 and form hydrogen bonding interactions with Ser190 and Gly214, (2) the phenyl of inhibitors can produce favorable van der Waals interactions with the residues His58, Cys191, Gln192, Trp211, Gly212, and Cys215. This systematic and comparative study can provide guidance for the choice of QM/MM MD methods and the designs of new potent inhibitors targeting trypsin.     

A computational study of the protein-ligand interactions in CDK2 inhibitors: using quantum mechanics/molecular mechanics interaction energy as a predictor of the biological activity

We report a combined quantum mechanics/molecular mechanics (QM/MM) study to determine the protein-ligand interaction energy between CDK2 (cyclin-dependent kinase 2) and five inhibitors with the N(2)-substituted 6-cyclohexyl-methoxy-purine scaffold. The computational results in this work show that the QM/MM interaction energy is strongly correlated to the biological activity and can be used as a predictor, at least within a family of substrates. A detailed analysis of the protein-ligand structures obtained from molecular dynamics simulations shows specific interactions within the active site that, in some cases, have not been reported before to our knowledge. The computed interaction energy gauges the strength of protein-ligand interactions. Finally, energy decomposition and multiple regression analyses were performed to check the contribution of the electrostatic and van der Waals energies to the total interaction energy and to show the capabilities of the computational model to identify new potent inhibitors.     

Recent developments of the quantum chemical cluster approach for modeling enzyme reactions

The quantum chemical cluster approach for modeling enzyme reactions is reviewed. Recent applications have used cluster models much larger than before which have given new modeling insights. One important and rather surprising feature is the fast convergence with cluster size of the energetics of the reactions. Even for reactions with significant charge separation it has in some cases been possible to obtain full convergence in the sense that dielectric cavity effects from outside the cluster do not contribute to any significant extent. Direct comparisons between quantum mechanics (QM)-only and QM/molecular mechanics (MM) calculations for quite large clusters in a case where the results differ significantly have shown that care has to be taken when using the QM/MM approach where there is strong charge polarization. Insights from the methods used, generally hybrid density functional methods, have also led to possibilities to give reasonable error limits for the results. Examples are finally given from the most extensive study using the cluster model, the one of oxygen formation at the oxygen-evolving complex in photosystem II.

QM/MM free energy simulations: recent progress and challenges

Due to the higher computational cost relative to pure molecular mechanical (MM) simulations, hybrid quantum mechanical/molecular mechanical (QM/MM) free energy simulations particularly require a careful consideration of balancing computational cost and accuracy. Here, we review several recent developments in free energy methods most relevant to QM/MM simulations and discuss several topics motivated by these developments using simple but informative examples that involve processes in water. For chemical reactions, we highlight the value of invoking enhanced sampling technique (e.g. replica-exchange) in umbrella sampling calculations and the value of including collective environmental variables (e.g. hydration level) in metadynamics simulations; we also illustrate the sensitivity of string calculations, especially free energy along the path, to various parameters in the computation. Alchemical free energy simulations with a specific thermodynamic cycle are used to probe the effect of including the first solvation shell into the QM region when computing solvation free energies. For cases where high-level QM/MM potential functions are needed, we analyse two different approaches: the QM/MM–MFEP method of Yang and co-workers and perturbative correction to low-level QM/MM free energy results. For the examples analysed here, both approaches seem productive although care needs to be exercised when analysing the perturbative corrections.     

Computer simulation based selection of optimal monomer for imprinting of tri-O-acetiladenosine in polymer matrix: vacuum calculations

Molecularly imprinted polymers can be anticipated as synthetic imitation of natural antibodies, receptors and enzymes. In case of successful imprinting the selectivity and affinity of the imprint for substrate molecules are comparable with those of natural counterparts. The selection of the optimal functional monomer, monomer/template ratio as well as choosing of polymerization solvent is crucial determinants of the successful imprinting. In the present study the simulation approach to the development of molecular imprinting polymers for the extraction of new protein kinase ATP-competitive inhibitors is presented. By imprinting tri-O-acetyladenosine into polymer matrix the synthetic reproduction of adenosine triphosphate binding site to protein kinases can be fabricated and further used for adenosine triphosphate analogs screening in different sources. The optimized geometrical structure and energy of the pre-polymerization complexes of tri-O-acetyladenosine (template) with three different monomers—methacrylic acid, 3-vinyl benzoic acid and acrylamide in vacuum were calculated using hybrid quantum mechanical/molecular mechanical (QM/MM) approach. These calculations demonstrate that methacrylic acid forms the most stable complex with template, the next is 3-vinyl benzoic acid complex and the third—acrylamide one. The bond energies of the complexes are shown to increase monotonically as more monomers are linked to the template. The same conclusions are made from purely quantum self-consistent field calculations of pre-polymerization complex energy and structure. Hybrid calculation is shown to be effective and can substantially accelerate the development of the imprinting technology.

Open image in new window

Figure Pre-polymerization complex of MIP with tri-O-acetiladenosine template with 5 metacrylic acid monomers

     

A QM/MM analysis of the conformations of crystalline sucrose moieties

Both ab initio quantum mechanics (QM) and molecular mechanics (MM) were used to produce a hybrid energy surface for sucrose that simultaneously provides low energies for conformations that are observed in crystal structures and high energies for most unobserved structures. HF/6-31G* QM energies were calculated for an analogue based on tetrahydropyran (THP) and tetrahydrofuran (THF). Remaining contributions to the potential energy of sucrose were calculated with MM. To do this, the MM surface for the analogue was subtracted from the MM surface for the disaccharide, and the QM surface for the analogue was added. Prediction of the distribution of observable geometries was enhanced by reducing the strength of the hydrogen bonding. Reduced hydrogen-bonding strength is probably useful because many crystalline sucrose moieties do not have intramolecular hydrogen bonds between the fructose and glucose residues. Therefore, hydrogen bonding does not play a large role in determining the molecular conformation. On the hybrid energy surface that was constructed with a dielectric constant of 3.5, the average potential energy of 23 sucrose moieties from crystal structures is 1.16 kcal/mol, and the population of observed structures drops off exponentially as the energy increases.