Mapping Interaction Energies in Chorismate Mutase with the Fragment

Subscriber access provided by University of Newcastle, Australia

Article

Mapping Interaction Energies in Chorismate Mutase with the Fragment Molecular Orbital Method Spencer R. Pruitt, and Casper Steinmann J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.6b12830 • Publication Date (Web): 08 Feb 2017 Downloaded from http://pubs.acs.org on February 12, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Mapping Interaction Energies in Chorismate Mutase with the Fragment Molecular Orbital Method Spencer R. Pruitt∗,† and Casper Steinmann∗,‡ †Academic & Research Computing, Worcester Polytechnic Institute, Worcester, MA 01602, United States ‡Centre for Computational Chemistry, School of Chemistry, University of Bristol, Bristol BS8 1TS, United Kingdom E-mail: [email protected]; [email protected]

1

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract The Claisen rearrangement of chorismate to prephenate is mapped across the entire reaction pathway using the fragment molecular orbital (FMO) method. Three basis sets (6-31G(d), cc-pVDZ, and pcseg-1) are studied to provide guidance towards obtaining high accuracy with the FMO method on such systems. Using a fragmentation scheme of one residue per fragment, the FMO method using the 6-31G(d) basis set and second-order Møller-Plesset perturbation theory (MP2) with the hybrid orbital projection (HOP) fragmentation scheme provides the most reliable results across the entire reaction pathway. Calculations using the multilayer FMO (ML-FMO) method are performed, and shown to be in agreement with single-layer calculations in all cases with differences of less than 2 kcal mol−1 for all tested basis set combinations along the entire reaction path. The use of restricted HartreeFock for the lower-level layer and MP2 for the higher-level layer gives the most consistent results when using the same basis set for both layers. Pair interaction energy decomposition analysis (PIEDA) calculations confirm that electrostatic interactions are the predominant force between three key arginine residues and chorismate, and that dispersion and charge transfer interactions in the binding pocket also play a role in the local chemistry of the reaction.

1. Introduction The study of chemical reactions in the liquid or solid state using traditional quantum mechanical (QM) methods is often unattainable due to the large scaling factor that accompanies the straightforward application of these methods. Several different strategies have been developed over the years to overcome this limitation, and have resulted in a wide variety of methodologies including hybrid quantum mechanics/molecular mechanics (QM/MM) 1,2 and numerous fragment based methods. 3,4 One such fragmentation method, the fragment molecular orbital (FMO) method, 5–7 has been extensively developed to treat a wide range 2


Page 2 of 42

Page 3 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of chemical and biological systems using all major QM levels of theory. Within the FMO formalism, a large system is fragmented into N non-overlapping fragments called monomers. Each fragment is treated using one of many available QM methods, with the interactions between the QM fragment and all other fragments represented through an electrostatic potential (ESP). The resulting self-consistent charge (SCC) iteration (i.e. the mutual convergence of all fragment densities) leads to a full N -body description of the electrostatic interactions. To capture short-range effects such as exchange-repulsion and charge-transfer, the QM properties of all fragment pairs (dimers) and optionally triples (trimers) are computed in the presences of the ESP of all remaining fragments. Although tremendous effort 8–13 has gone into making the FMO method applicable to optimizing geometries 14–16 and transition states, 17,18 the FMO method has historically been applied as an energy-refinement method 7,19,20 in ligand binding studies, 21–26 reaction barrier estimation, 27 or identification of key residues through inspection of the interaction energies between substrates and surrounding residues. 28,29 Although the FMO method is significantly more capable and efficient than performing a full QM calculation on a large molecular system, 4,30 accurate mapping of catalytic reactions in many cases requires the inclusion of the environment (solvent, ions, etc). To provide an even more efficient method to map out complete reactions, the effective fragment molecular orbital 31–33 (EFMO) method was devised. The EFMO method merges the effective fragment potential (EFP) method 34,35 and the FMO method, lowering the computational cost even further by computing fragments in the gas phase, and including polarization effects with the ab initio based EFP method. The approximate gradients initially developed for the EFMO method enabled the mapping of the Claisen rearrangement of chorismate to prephenate 36–38 which is catalyzed by Bacillus Subtilis chorismate mutase. 39 The chorismate to prephenate reaction (Scheme 1) has an experimental activation free energy of ∆G‡ = 15.4 kcal mol−1 in the enzyme, and a barrier of ∆G‡ = 24.5 kcal mol−1 in aqueous solution. 40 In both cases the reaction is believed to be a unimolecular pericyclic reaction. The reaction has garnered much attention from

3



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 42

both experimental and computational chemists to explain the origin of the catalytic effect in the enzyme. Mutation studies have identified Arg90 as playing a major role in transition O-

O-

O

O O O-

4

O

3

3

-

O-

2O 1 O

3 4

4

O

O

2

1

2

O

O

1

O-

OH

HO OH

1

2

3

Scheme 1: Illustration of the conversion of chorismate (1) through the transition state (2) to prephenate (3). Figure is reproduced from Ref 37 with permission. state stabilization during the reaction. 41,42 The reaction has also been investigated in detail with semi-empirical QM/MM methodologies, and subsequently refined using high-level correlated methods 43,44 which also verified the stabilizing effect of the enzyme. 45 In this work we hope to shed additional light on important interactions between nearby residues and the substrate along the entire reaction pathway using the FMO method. FMO offers a unique capability in that it can provide insight into energy contributions from each fragment in ligand-binding studies. Besides pair-interaction energies (PIEs) between fragments, the interaction energy can be decomposed into four terms using the pair interaction energy decomposition analysis (PIEDA). 46 The ability to perform an energy decomposition analysis 47,48 (EDA) in drug-discovery applications 49 has become vital, as was recently highlighted by Heifetz et al. 50 The present study for the first time applies the PIEDA procedure along an entire reaction coordinate. As chorismate mutase converts chorismate to prephenate different energy contributions to the total interaction energy are investigated, and insights from such a detailed analysis may prove relevant to the future design of new drug-like molecules, or investigation of different catalytic pathways. To provide guidance for the use of PIEDA during such investigations, a variety of possible QM method and basis set combinations are explored, with the goal of providing a benchmark for computation of reaction barriers that provides both low 4


Page 5 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


cost and high accuracy. Calculated enthalpy barriers and averaged PIEs, estimated using multiple reaction trajectories, are presented. The paper is organized as follows: First, the theoretical background for the FMO method (including PIEDA) is outlined, followed by the computational details. Second, results concerning energy barriers at different levels of theory and approximations within the FMO method are presented. Next, a detailed analysis of the PIEDA results is presented for all tested reaction trajectories. The paper concludes with a summary of the findings, as well as the current state of and growing potential for use of the FMO method in drug discovery research.

2. Theoretical Methods 2.1 The Fragment Molecular Orbital Method The complete formalism of the FMO method, including rigorous derivations, has been described previously elsewhere. 3 The formalism described herein is therefore only a brief outline of the equations relevant to the present study. The total energy of a system of N fragments using the two-body FMO method (FMO2) can be written as 51 E

FMO2

=

N X

EI′

I

+

N X

∆EIJ

(1)

I>J

where the PIEs (second term in eq 1) are defined as

′ ∆EIJ = EIJ − EI′ − EJ′ + Tr ∆DIJ VIJ

(2)

′ and the energies EX are internal energies of monomers (I) or dimers (IJ)

′ EX = EX − Tr DX VX .

5


(3)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 42

The addition of fragment triples (IJK) to the energy expression (EIJK ) would give the total three-body FMO method (FMO3) energy. 52,53 The dimer density difference matrix ∆DIJ , defined as

∆DIJ = DIJ − (DI ⊕ DJ ),

(4)

represents the density difference between a pair of fragments (IJ) and the direct sum of the individual fragment densities of I and J. VX is the ESP ”felt” by a monomer (I) or a dimer (IJ) from all surrounding fragments. While the ESP (VX ) can be calculated exactly 5,6 from fragment densities and two-electron integrals, the computation requires a non-negligible computational cost on the order of O(N 2 ). To reduce this cost, the energy expression presented in eq 1 was formulated to allow for approximations to the ESP based on Mulliken charges 51 or potential derived charges. 54 In this work no attempt is made to interfere with the delicate balance of the ESP. However, it should be noted that the ESP can have a very strong influence on the accuracy of the energies obtained, as discussed in detail in Ref 54. For bonded molecular systems treated within the FMO formalism, covalent bonds between fragments are detached in such a way that one fragment is assigned the entire bond, and the other fragment is assigned no part of the bond. After detachment, the boundaries between fragments are handled using either the hybrid orbital projection (HOP) method 6 or the adaptive frozen orbital 55 (AFO) method. Both of these bond fractionation techniques seek to contain the electronic density within a given fragment. Augmenting the FMO2 reference energy with correlation energy calculated using secondorder Møller-Plesset perturbation theory (MP2) 56,57 is achieved in a similar fashion to the FMO2 energy E

FMO2−MP2

=E

FMO2

+

N X I

6

EIcorr

+

N X I>J


corr ∆EIJ

(5)

Page 7 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


where corr corr ∆EIJ = EIJ − EIcorr − EJcorr

(6)

int The total interaction energy within a dimer (IJ), ∆EIJ , can be divided into reference and

correlated contributions as, int corr ∆EIJ = ∆EIJ + ∆EIJ .

(7)

The correlation energy can also be evaluated at the three-body level (FMO3-MP2). 57

2.2 Approximations The above formulation of the FMO2-MP2 energy allows for nearly linear scaling with respect to system size, provided that dimer energies are approximated at long distances. 58,59 In addition to explicit distance based approximations to the fragment interactions, additional reductions in computational cost can be achieved through approximations to the ESP, 51 or through the use of the multilayer formulation of FMO 60 (ML-FMO). Within the ML-FMO formalism, each fragment is assigned to a layer, with each layer described using a different basis set and/or QM level of theory. Inter-layer boundaries are straight forwardly treated, as they coincide with inter-fragment boundaries, with dimers only evaluated if the fragments that compose the dimer are included in the same layer. All fragments in a higher-level layer are also included in all lower-level layers. The standard application of the ML-FMO method is to use a lower-level of theory to describe a surrounding region, and embedding a higher-level region within the ESP generated from the lower-level of theory.

2.3 Pair Interaction Energy Decomposition Analysis The interaction energy within a dimer (IJ) in (eq 7) can be decomposed as 46,49 disp CT+mix int es ex . + ∆EIJ ∆EIJ = ∆EIJ + ∆EIJ + ∆EIJ

7


(8)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 42

where terms 1-4 on the rhs of eq 8 represent electrostatic, exchange repulsion, charge transfer, and dispersion interactions between fragments in a given dimer IJ. For systems with covalent bonds, the expression in eq 8 includes additional terms. 46 However, the present study is focused on the interaction energy between chorismate and the surrounding protein, which share no covalent bonds, requiring that no additional terms be considered. The focus of the present study is to investigate the accuracy of various basis sets across an entire reaction pathway, as well as determine what insights into the chorismate reaction can be obtained from PIEDA calculations when applied to an entire reaction path. In order to obtain the highest accuracy from FMO calculations on peptide based systems in terms of absolute energies (i.e. when not performing PIEDA calculations) one of the following fragmentation schemes should be employed: two residues per fragment using FMO2, or one residue per fragment using FMO3. 54 Based on the outlined focus of this study, a fragmentation scheme of one residue per fragment is used in all calculations.

3. Computational Details All calculations in this study were performed using the FMO method as implemented in the GAMESS 61 computational chemistry software package. For each type of calculation performed, the following basis sets were used: STO-3G, 62 6-31G(d), 62,63 cc-pVDZ, 64,65 and the polarization consistent segment contracted basis sets (pcseg-n). 66,67 All structures were obtained from previous work by Steinmann et al 37 who calculated reaction barriers for the conversion of chorismate to prephenate in chorismate mutase obtained at the EFMO 31,36 RHF/6-31G(d) level of theory. The reaction coordinate R is defined as the difference in the bond length between the breaking O2-C1 bond and the forming C4-C3 bond (see atom numbers in Scheme 1): R = R21 − R43 .

8


(9)

Page 9 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


This definition has a reaction coordinate of R = −2.0 ˚ A for chorismate (R21 = 1.4 ˚ A and R43 = 3.4 ˚ A) and correspondingly for prephenate at R = 2.1 ˚ A (R21 = 3.3 ˚ A and R43 = 1.4 ˚ A). All structures were prepared using FragIt, 68 which has been updated to include support for the HOP fragmentation scheme. Default values were used for all distance based cutoffs in FMO, including approximations to the ESP, separated dimers, and neglect of correlation in well-separated dimers. Hybrid molecular orbitals for the STO-3G, 6-31G(d), and cc-pVDZ basis sets used in FMO calculations involving the HOP fragmentation scheme are freely distributed with GAMESS. However, for the pcseg-n basis sets HMOs were constructed using the following procedure: 1. The geometry of methane was optimized using C3v symmetry with one C – H bond directed along the principle axis (z-axis in GAMESS) at the RHF level of theory. 2. After geometry optimization, the resulting molecular orbitals at the stationary point were localized using the Edminston-Ruedenberg localization 69 scheme. 3. The localized orbitals for the carbon atom were identified and extracted for use in the computations at the fragmentation points. The resulting HMOs used in the input files can be found in the supporting information.

3.1 Models A set of representative models for use in FMO calculations were prepared. Reaction barriers were calculated at the MP2 level of theory, as previous studies have shown that some level of electron correlation is required. 37,44 All models used in the calculation of reaction barriers are presented in Table 1. Standard FMO2 calculations with MP2 (FMO2-X/B) are designated Rn, with energies calculated according to eq 5 with electronic structure method X and basis set B. For example, R1 is a FMO2 calculation using MP2 to calculate energies for both monomers and dimers using the 6-31G(d) basis set, i.e. FMO2-MP2/6-31G(d). 9



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment

Page 10 of 42

Page 11 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


layer (all residues and water molecules with at least one atom within 4.5 ˚ A of chorismate). The remainder of the protein, shown in gray, is treated as the lower-level layer. Table 1: Different levels of theory used for the FMO models in the calculation of reaction barriers. Regular FMO methods are designated as models Rn, and ML-FMO methods are designated as models Mn. Model R1 R2 R3 R4 R5 R6 M1 M2 M3 M4 M5 M6 M7 M8 M9 M10 M11 M12 M13 M14 M15 M16

Method FMO2-MP2 FMO2-MP2 FMO2-MP2 FMO2-MP2 FMO2-MP2 FMO2-MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2 FMO2-RHF:MP2

Basis 6-31G(d) 6-31G(d) pcseg-1 pcseg-1 cc-pVDZ cc-pVDZ 6-31G(d):6-31G(d) 6-31G(d):6-31G(d) STO-3G:6-31G(d) STO-3G:6-31G(d) STO-3G:cc-pVDZ STO-3G:cc-pVDZ 6-31G(d):cc-pVDZ 6-31G(d):cc-pVDZ pcseg-0:6-31G(d) pcseg-0:6-31G(d) pcseg-0:pcseg-1 pcseg-0:pcseg-1 pcseg-1:pcseg-1 pcseg-1:pcseg-1 cc-pVDZ:cc-pVDZ cc-pVDZ:cc-pVDZ

11


Bonds HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO HOP AFO



Page 12 of 42

Page 13 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the FMO2-MP2/AFO results. The qualitative difference between using the HOP or AFO schemes is not large when inspecting Figures 2A and B. However, Figure 2C shows the difference between the reaction barriers (E(Rn) − E(Rm)) computed using the HOP scheme and those computed using the AFO scheme increases as the reaction progresses towards the transition state (R = −0.57 ˚ A) . The maximum difference for all tested double zeta quality basis sets at R = −0.3 ˚ A occurs using model R3, with a calculated difference of 2.5 kcal mol−1 . As the reaction progresses from the transition state towards the product, the differences between the HOP and AFO schemes flatten out with the R3 model becoming constant at around 1.5 kcal mol−1 . Models R1 and R5 have comparable deviations near the transition state of approximately 1 kcal mol−1 lower than R3. After progressing past the transition state, the differences between the HOP and AFO schemes for both R1 and R5 flatten out and becomes constant at around 1 kcal mol−1 in the same manner as R3. Table 2: Reaction barrier heights in kcal mol−1 for FMO2-MP2 and FMO3-MP2.

6-31G(d) cc-pVDZ pcseg-1

HOP FMO2 FMO3 16.4 16.8 11.5 16.1 12.8 16.2

AFO FMO2 FMO3 14.8 17.0 10.0 16.4 10.5 16.6

As was noted in section 2.2, high-accuracy results can be obtained by using either FMO3 with one residue per fragment, or FMO2 with two residues per fragment. In Table 2 results are shown for computed barrier heights at the FMO2-MP2 and FMO3-MP2 levels of theory for comparison. For FMO3-MP2/6-31G(d) the computed barrier heights are 16.8 and 17.0 kcal mol−1 for the HOP and AFO based calculations, respectively, which shows excellent agreement when compared to FMO2-MP2 for the HOP scheme with a deviation of only 0.4 kcal mol−1 . However, calculations using the AFO scheme deviate by 2.2 kcal mol−1 . The larger basis sets shows deviations of 4.6 and 3.4 kcal mol−1 for cc-pVDZ and pcseg-1, respectively for the HOP results. In summary, the FMO2-MP2/6-31G(d) calculations compare favorably with the FMO313



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

MP2 calculations using the HOP scheme and any of the three basis sets. In light of recommendations regarding FMO accuracy as outlined above, it is important to remember the effect of the ESP on the accuracy. Previous studies have shown that small imbalances between the ESP and the level of many-body decomposition (e.g. FMO2 or FMO3) can cause large reductions in accuracy of the overall method. 70 For example, the use of a large basis set (e.g. cc-pVDZ) with the FMO2 method can reduce the accuracy such that the use of the FMO3 method containing higher order many-body interactions is required to obtain the same accuracy as smaller basis set (e.g. 6-31G(d)) FMO2 calculations.

4.2 Multilayer Calculations In this section the ability of the ML-FMO method to recover the computed reaction barriers calculated using the standard FMO2-MP2 method is investigated. When employing the MLFMO method, a region of the system of interest is acknowledged to be of greater importance than the rest of the system, with only this important region treated at the designated higherlevel of theory. The rest of the system is treated at the designated lower-level of theory in an effort to reduce computational time without compromising accuracy (see Figure 1). Figure 3 shows the difference between the reference calculations and the ML calculations (E(Rn) − E(Mm)) as a function of the reaction coordinate. Figure 3A shows results at the FMO2-MP2/6-31G(d) level of theory (models R1 and R2 for HOP and AFO, respectively) as the reference values to quantify the performance of ML-FMO calculations that use 631G(d) (models M1 and M2), STO-3G (models M3 and M4), and pcseg-0 (models M9 and M10) basis sets in the lower-level layer, and the 6-31G(d) basis set in the higher-level layer. In Figure 3B, FMO2-MP2/cc-pVDZ (models R5 and R6) are used as the reference calculations to quantify how well FMO2-RHF:MP2 ML calculations perform that use one of the STO-3G (models M5 and M6), 6-31G(d) (models M7 and M8), or cc-pVDZ (models M15 and M16) basis sets in the lower-level layer, and the cc-pVDZ basis set in the higherlevel layer. And in Figure 3C, the accuracy of using either the pcseg-0 or pcseg-1 basis set in 14


Page 14 of 42

Page 15 of 42




1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Model M4 (dashed orange) follows the behavior of M3 very closely. As was observed for the results in Figure 3A, Figure 3B shows that the combination with the cc-pVDZ:cc-pVDZ basis set combination set in both the higher- and lower-level layers (models M15 and M16) gives results that deviate only by 0.34 kcal mol−1 around the transition state when comparing with the corresponding FMO2-MP2/cc-pVDZ reaction barrier. For example, when using the AFO scheme (dashed green) the reference energy gives E TS (M7) = 11.54 kcal mol−1 , whereas the similar multilayer calculation gives E TS (M15) = 11.30 kcal mol−1 . The same small deviation can be observed around the transition state when using the cc-pVDZ:cc-pVDZ basis set combination with the HOP scheme (solid green), mirroring the behavior of the 6-31G(d):6-31G(d) results presented in Figure 3A. A reduction in the size of the basis set used to describe the lower-level layer shows a dramatic change in the accuracy of the ML-FMO method when using the cc-pVDZ basis set for the higherlevel layer. Models M7 and M8 over estimate the region of the reaction pathway around the transition state by over 1 kcal mol−1 for both the HOP based (solid blue) and AFO based (dashed blue) calculations. In a similar fashion, the use of the STO-3G basis set overestimates the same region around the transition state by 1.5 kcal/mol. Models M5 and M6 produce nearly identical results, yet also fail to reproduce the latter half of the reaction when compared to the reference values for both HOP and AFO based calculations respectively. When comparing pcseg-0:pcseg-1 calculations with FMO2-MP2/pcseg-1 (R5 and R6) an increased deviation is observed, but contrary to the results in Figure 3A and B the deviation to the reference barrier is in general much less pronounced when using the segment contracted polarization consistent basis sets. A 0.4 kcal mol−1 overestimation of the barrier height can also be observed using pcseg-0:pcseg-1, which is consistent between both the HOP and AFO schemes. As the reaction moves towards the product, the reaction energy is underestimated by 0.8 kcal mol−1 for HOP and 1.0 kcal mol−1 for AFO. This trend is similar to, but less pronounced than, the results for 6-31G(d):cc-pVDZ in Figure 3B. The application of the

16


Page 16 of 42

Page 17 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


pcseg-1 basis set for the lower- and higher-level layers (models M13 and M14) produces more consistent errors across the entire reaction pathway. A deviation around the transition state of 0.5 kcal mol−1 is present using the HOP scheme (solid black). However, the same variation is absent for the corresponding AFO calculations. Comparing all of the ML calculations performed, none of the models deviate by more than 2 kcal mol−1 from the reference calculations for all tested basis sets and fragmentation schemes. However, not all combinations produced consistent errors across the entire reaction pathway, which is vital for calculation of accurate reaction enthalpies and barriers. As shown in Figures 3A and 3B, the use of a basis set for the lower-level layer that is too small in comparison to the higher-level layer basis set can cause variability in the accuracy of the ML-FMO method at different points along the reaction pathway. Calculations based on the ML-FMO method that use the same basis set for the lower- and higher-level layers (models M1, M2, M15, and M16) produce consistent results across the entire reaction pathway, and are a viable option to reduce computational cost through the reduction of the QM level of theory used on the lower-level layer of the system of interest. Results in Figure 3C show that the ML-FMO method produces more consistent results with the use of the smaller the segment contracted polarization consistent basis set for the lower-level layer (models M11 and M12), in contrast to the use of either STO-3G or 6-31G(d) with a larger basis set for the higher-level layer (models M3 through M8). The results produced by the segment contracted polarization consistent basis sets compare well with the much more expensive cc-pVDZ:cc-pVDZ calculations in terms of reproducing FMO2-MP2 results at either the cc-pVDZ or pcseg-1 level of theory. Even when using the segment contracted polarization consistent basis sets, the use of the same basis set for both lower- and higher-level layers (models M13 and M14) in ML-FMO calculations still produce more consistent results across the entire reaction pathway.

17




Page 18 of 42

Page 19 of 42




1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

that the electrostatic interaction between chorismate and Arg90 is crucial for transition state stabilization (15 kcal mol−1 ), in agreement with previous studies. 71,72 Arg7 also contributes 7 kcal mol−1 to the stabilization of the transition state, while the contribution from Arg63 is less than half that of Arg7, approximately 3 kcal mol−1 at the transition state geometry. As the reaction continues towards the product, the electrostatic interaction between Arg90 and chorismate returns to the same relative magnitude as the initial reactant interactions. In contrast, the interaction between Arg7 and chorismate decreases significantly, dropping by nearly 10 kcal mol−1 . While Arg63 does not appear to play an important role in the electrostatic interaction between Arg63 and chorismate at the transition state, there is a net stabilization of nearly 10 kcal mol−1 between Arg63 and chorismate as the reaction progresses towards the product. Upon inspection of key distances (Figure 4), it can be observed that the main hydrogen bond between chorismate and Arg63 undergoes a contraction from 1.8 ˚ A to 1.6 ˚ A over the course of the reaction. The role that dispersion interactions play is shown in Figure 5B. The general trend observed for the electrostatic interactions between the three arginine residues and chorismate is also apparent for the dispersion interactions. Dispersion also contributes to transition state stabilization, with the strongest net stabilization of the interaction energy of approximately 1.5 kcal mol−1 at the transition state arising from the Arg90/chorismate interaction. However while Arg7 does contribute to stabilization of the transition state electrostatic energy, Arg7 does not contribute to transition state stabilization in terms of dispersion effects. The interaction between Arg7 and chorismate does show a net reduction in dispersion interaction of around 0.5 kcal mol−1 , which mirrors the behavior of the electrostatic interaction between Arg7 and chorismate. Arg63 shows a slight increase in the dispersion interaction with chorismate across the entire reaction of approximately 1 kcal mol−1 . Figure 5C shows the relative charge-transfer interaction energies across the entire reaction pathway. As was observed for both electrostatic and dispersion interactions, Arg90 plays the most important role in terms of stabilization of the charge transfer energy around the

20


Page 20 of 42

Page 21 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


transition state with a net increase of approximately 1.3 kcal mol−1 . The behavior of the charge transfer energy between Arg7 and chorismate is markedly different than the behavior observed for electrostatics and dispersion. In both Figure 5A and Figure 5B, the Arg7 residue provides a stabilization of the transition state at nearly the same minimum on the reaction coordinate. However in Figure 5C, the charge transfer interaction shows a clear stabilizing effect after the observed transition state. This shift in stabilization could aid in shifting the reaction past the transition state and towards the products. Combine with the rapid increase in charge transfer stabilization between Arg63 and chorismate, this comparatively small interplay between chorismate and the Arg7 and Arg63 residues appears to play an important role in the chemistry of the reaction pathway. Finally, the relative exchange-repulsion interaction energies between the three arginines are shown in Figure 5D. Similarly to the other panels in Figure 5, Arg63 is observed to behave differently. Indeed the exchange repulsion term is increased at the product state by almost 10 kcal mol−1 compared to the reactant state. This further confirms that the distance between Arg63 and chorismate is decreased since the exchange-repulsion term is short-ranged. On the other hand, both Arg7 and Arg90 exhibit similar drops in the exchangerepulsion energy confirming that the distance between those two residues and chorismate is increased as the reaction goes from reactant to product. 4.3.2 Total Interaction Energies The investigation of the interaction energies concludes with a discussion of the net change in the total decomposed interaction between the reactant and the product. Figure 6 shows the net change of each energy component for all trajectories included in this study, and an averaged result over all trajectories. Overall, electrostatic interactions at the product state decrease by an average of 15 kcal mol−1 compared to the reactant state, but the variation between trajectories is as low as 10 kcal mol−1 for trajectory 7, and upwards of 21 kcal mol−1 for trajectory 5. It can also be observed that when the net change in the electrostatic

21




Page 22 of 42

Page 23 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the FMO method for mapping reaction pathways in drug discovery studies. In addition, the PIEDA as part of the FMO method was studied as an efficient mechanism for obtaining insights into the local chemistry of such reaction pathways in biological systems. Three basis sets (6-31G(d), cc-pVDZ, and pcseg-1) were studied in conjunction with the ML-FMO method. While no qualitative differences between the two fragmentation schemes available in the FMO method (HOP or AFO) were observed, quantitative differences were shown to exist in the computed barriers that are dependent on the basis set. To compute accurate reaction energy barriers using the FMO2-MP2 method, it was found that the 631G(d) basis set with the HOP fragmentation scheme was the most accurate compared to FMO3-MP2 results. Although the AFO scheme is slightly easier to use for the end-user, the computed energies using the AFO scheme underestimate the reaction barriers using the HOP scheme at the FMO2 level of theory, possibly due to the way polarization is restricted across covalent bonds. To help end-users make use of the HOP scheme, support for the use of the hybrid orbitals was added to the FragIt program. It was determined that while the use of large basis sets in quantum chemistry is viewed as a way to provide greater accuracy, within the FMO formalism the balance of the environmental ESP with the fully QM fragment calculations must be taken into account. Using a fragmentation scheme of one residue per fragment, the 6-31G(d) basis set at the FMO2MP2 level of theory provides the most reliable results across the entire reaction pathway. The accuracy obtained with modern double-zeta basis sets (e.g. cc-pVDZ) requires either FMO3-MP2, or adjusted settings for the ESP to be competitive with FMO2-MP2/6-31G(d). It is possible, but not studied here, that the role of the embedding ESP could change the accuracy of the results presented, and is worthy of an independent study. The FMO method is a delicate balance of approximations and error-cancellation, and the ML-FMO method is no exception. The employed approximations must co-exist with the error-cancellation in order to obtain accurate results. Calculations using the ML-FMO method were performed, and the agreement with single-layer calculations were in all cases

23



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 42

below 2 kcal mol−1 for all tested basis set combinations. It was found that the use of RHF for the lower-level layer and MP2 for the higher-level layer gives the most consistent results when using the same basis set for both layers. The reduction in computational cost can be beneficial when performing explorative calculations on potential energy surfaces. PIEDA calculations confirm that all three arginine residues interact with chorismate predominantly through an electrostatic interaction (Figure 5A). Despite the large contribution from electrostatics to the overall reaction, the comparatively small energetics of dispersion and charge transfer between chorismate and the Arg7 and Arg63 residues also play an important role in the local chemistry of the reaction. While the FMO method is still far from what can be considered a ”black box” method, the development of tools such as FragIt, and the ongoing addition of functionality, has allowed the FMO method to potentially be used by non-experts in the field. Continued use of fragmentation methods, such as the FMO method, will add to the growing body of literature, giving guidance to the proper use, failings, and best practices for this unique class of methods.

Acknowledgement C.S. thanks the Danish Council for Independent Research (the Sapere Aude program) for financial support (Grant No. 4181-00370) and Dr. Kara Ranaghan for fruitful discussions.

Supporting Information Available A document and a zip-archive both containing the following files: Hybrid orbital projector orbitals used in FMO calculations. FMO input and output files are available on Figshare (url: https://figshare.com/projects/Mapping_Interaction_Energies_in_Chorismate_Mutase_ with_the_Fragment_Molecular_Orbital_Method/18580) This material is available free of charge via the Internet at http://pubs.acs.org/. 24


Page 25 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


References (1) Senn, H. M.; Thiel, W. QM/MM Methods for Biomolecular Systems. Angew. Chem. Int. Ed. 2009, 48, 1198–1229. (2) van der Kamp, M. W.; Mulholland, A. J. Combined Quantum Mechanics/Molecular Mechanics (QM/MM) Methods in Computational Enzymology. Biochemistry 2013, 52, 2708–2728. (3) Gordon, M. S.; Fedorov, D. G.; Pruitt, S. R.; Slipchenko, L. V. Fragmentation Methods: A Route to Accurate Calculations on Large Systems. Chem. Rev. 2012, 112, 632–672. (4) Pruitt, S. R.; Bertoni, C.; Brorsen, K. R.; Gordon, M. S. Efficient and Accurate Fragmentation Methods. Accounts Chem. Res. 2014, 47, 2786–2794. (5) Kitaura, K.; Ikeo, E.; Asada, T.; Nakano, T.; Uebayasi, M. Fragment Molecular Orbital Method: An Approximate Computational Method for Large Molecules. Chem. Phys. Lett. 1999, 313, 701–706. (6) Nakano, T.; Kaminuma, T.; Sato, T.; Akiyama, Y.; Uebayasi, M.; Kitaura, K. Fragment Molecular Orbital Method: Application to Polypeptides. Chem. Phys. Lett. 2000, 318, 614–618. (7) Fedorov, D. G.; Kitaura, K. Extending the Power of Quantum Chemistry to Large Systems with the Fragment Molecular Orbital Method. J. Phys. Chem. A 2007, 111, 6904–6914. (8) Nagata, T.; Brorsen, K.; Fedorov, D. G.; Kitaura, K.; Gordon, M. S. Fully Analytic Energy Gradient in the Fragment Molecular Orbital Method. J. Chem. Phys. 2011, 134, 124115. (9) Nagata, T.; Fedorov, D. G.; Kitaura, K. Analytic Gradient for the Embedding Potential

25



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

with Approximations in the Fragment Molecular Orbital Method. Chem. Phys. Lett. 2012, 544, 87–93. (10) Nagata, T.; Fedorov, D. G.; Li, H.; Kitaura, K. Analytic Gradient for Second Order Møller-Plesset Perturbation Theory with the Polarizable Continuum Model Based on the Fragment Molecular Orbital Method. J. Chem. Phys. 2012, 136, 204112. (11) Brorsen, K. R.; Minezawa, N.; Xu, F.; Windus, T. L.; Gordon, M. S. Fragment Molecular Orbital Molecular Dynamics with the Fully Analytic Energy Gradient. J. Chem. Theory Comput. 2012, 8, 5008–5012. (12) Brorsen, K. R.; Zahariev, F.; Nakata, H.; Fedorov, D. G.; Gordon, M. S. Analytic Gradient for Density Functional Theory Based on the Fragment Molecular Orbital Method. J. Chem. Theory Comput. 2014, 10, 5297–5307. (13) Nakata, H.; Fedorov, D. G.; Yokojima, S.; Kitaura, K.; Nakamura, S. Derivatives of the Approximated Electrostatic Potentials in Unrestricted Hartree–Fock Based on the Fragment Molecular Orbital Method and an Application to Polymer Radicals. Theor. Chem. Acc. 2014, 133 . (14) Fedorov, D. G.; Alexeev, Y.; Kitaura, K. Geometry Optimization of the Active Site of a Large System with the Fragment Molecular Orbital Method. J. Phys. Chem. Lett. 2011, 2, 282–288. (15) Nakata, H.; Fedorov, D. G.; Nagata, T.; Kitaura, K.; Nakamura, S. Simulations of Chemical Reactions with the Frozen Domain Formulation of the Fragment Molecular Orbital Method. J. Chem. Theory Comput. 2015, 11, 3053–3064. (16) Nakata, H.; Fedorov, D. G. Efficient Geometry Optimization of Large Molecular Systems in Solution Using the Fragment Molecular Orbital Method. J. Phys. Chem. A. 2016,

26


Page 26 of 42

Page 27 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(17) Nakata, H.; Fedorov, D. G.; Yokojima, S.; Kitaura, K.; Nakamura, S. Efficient Vibrational Analysis for Unrestricted Hartree–Fock Based on the Fragment Molecular Orbital Method. Chem. Phys. Lett. 2014, 603, 67–74. (18) Nakata, H.; Nagata, T.; Fedorov, D. G.; Yokojima, S.; Kitaura, K.; Nakamura, S. Analytic Second Derivatives of the Energy in the Fragment Molecular Orbital Method. J. Chem. Phys. 2013, 138, 164103. (19) Fedorov, D. G.; Nagata, T.; Kitaura, K. Exploring Chemistry with the Fragment Molecular Orbital Method. Phys. Chem. Chem. Phys. 2012, 14, 7562. (20) Fedorov, D. The Fragment Molecular Orbital Method: Practical Applications to Large Molecular Systems; CRC Press/Taylor & Francis: Boca Raton, 2009. (21) Poongavanam, V.; Steinmann, C.; Kongsted, J. Inhibitor Ranking through QM Based Chelation Calculations for Virtual Screening of HIV-1 RNase H Inhibition. PLoS ONE 2014, 9, e98659. (22) Jensen, J. H.; Willemo¨es, M.; Winther, J. R.; Vico, L. D. In Silico Prediction of Mutant HIV-1 Proteases Cleaving a Target Sequence. PLoS ONE 2014, 9, e95833. (23) Ito, M.; Brinck, T. Novel Approach for Identifying Key Residues in Enzymatic Reactions: Proton Abstraction in Ketosteroid Isomerase. J. Phys. Chem. B. 2014, 118, 13050–13058. (24) Ishikawa, T.; Burri, R. R.; Kamatari, Y. O.; Sakuraba, S.; Matubayasi, N.; Kitao, A.; Kuwata, K. A Theoretical Study of the Two Binding Modes Between Lysozyme and Tri-NAG with an Explicit Solvent Model Based on the Fragment Molecular Orbital Method. Phys. Chem. Chem. Phys. 2013, 15, 3646. (25) Yoshino, R.; Yasuo, N.; Inaoka, D. K.; Hagiwara, Y.; Ohno, K.; Orita, M.; Inoue, M.; Shiba, T.; Harada, S.; Honma, T. et al. Pharmacophore Modeling for Anti-Chagas 27



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Drug Design Using the Fragment Molecular Orbital Method. PLOS ONE 2015, 10, e0125829. ¨ (26) Oberg, H.; Brinck, T. Fragment Molecular Orbital Study of the cAMP-Dependent Protein Kinase Catalyzed Phosphoryl Transfer: A Comparison with the Differential Transition State Stabilization Method. Phys. Chem. Chem. Phys. 2016, 18, 15153– 15161. (27) Polyakov, I. V.; Grigorenko, B. L.; Moskovsky, A. A.; Pentkovski, V. M.; Nemukhin, A. V. Towards Quantum-Based Modeling of Enzymatic Reaction Pathways: Application to the Acetylholinesterase Catalysis. Chem. Phys. Lett. 2013, 556, 251– 255. (28) Nakamura, T.; Yamaguchi, A.; Kondo, H.; Watanabe, H.; Kurihara, T.; Esaki, N.; Hirono, S.; Tanaka, S. 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. J. Comput. Chem. 2009, 30, 2625–2634. (29) Hitaoka, S.; Chuman, H.; Yoshizawa, K. A QSAR Study on the Inhibition Mechanism of Matrix Metalloproteinase-12 by Arylsulfone Analogs Based on Molecular Orbital Calculations. Org. Biomol. Chem. 2015, 13, 793–806. (30) Gordon, M. S.; Mullin, J. M.; Pruitt, S. R.; Roskop, L. B.; Slipchenko, L. V.; Boatz, J. A. Accurate Methods for Large Molecular Systems. J. Phys. Chem. B 2009, 113, 9646–9663. (31) Steinmann, C.; Fedorov, D. G.; Jensen, J. H. Effective Fragment Molecular Orbital Method: A Merger of the Effective Fragment Potential and Fragment Molecular Orbital Methods †. J. Phys. Chem. A 2010, 114, 8705–8712. (32) Pruitt, S. R.; Steinmann, C.; Jensen, J. H.; Gordon, M. S. Fully Integrated Effective Fragment Molecular Orbital Method. J. Chem. Theory Comput. 2013, 9, 2235–2249. 28


Page 28 of 42

Page 29 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(33) Steinmann, C.; Jensen, J. H. In Fragmentation: Toward Accurate Calculations on Complex Molecular Systems, 1st ed.; Gordon, M. S., Ed.; John Wiley & Sons, Ltd., 2017; Chapter 5. (34) Day, P. N.; Jensen, J. H.; Gordon, M. S.; Webb, S. P.; Stevens, W. J.; Krauss, M.; Garmer, D.; Basch, H.; Cohen, D. An Effective Fragment Method for Modeling Solvent Effects in Quantum Mechanical Calculations. J. Chem. Phys. 1996, 105, 1968–1986. (35) Gordon, M. S.; Freitag, M. A.; Bandyopadhyay, P.; Jensen, J. H.; Kairys, V.; Stevens, W. J. The Effective Fragment Potential Method:

A QM-Based MM Ap-

proach to Modeling Environmental Effects in Chemistry. J. Phys. Chem. A 2001, 105, 293–307. (36) Steinmann, C.; Fedorov, D. G.; Jensen, J. H. The Effective Fragment Molecular Orbital Method for Fragments Connected by Covalent Bonds. PLoS ONE 2012, 7, e41117. (37) Steinmann, C.; Fedorov, D. G.; Jensen, J. H. Mapping Enzymatic Catalysis Using the Effective Fragment Molecular Orbital Method: Towards all ab initio Biochemistry. PLoS ONE 2013, 8, e60602. (38) Christensen, A. S.; Steinmann, C.; Fedorov, D. G.; Jensen, J. H. Hybrid RHF/MP2 Geometry Optimizations with the Effective Fragment Molecular Orbital Method. PLoS ONE 2014, 9, e88800. (39) Chook, Y. M.; Ke, H.; Lipscomb, W. N. Crystal Structures of the Monofunctional Chorismate Nutase from Bacillus Subtilis and its Complex with a Transition State Analog. Proc. Natl. Acad. Sci. 1993, 90, 8600–8603. (40) Kast, P.; Asif-Ullah, M.; Hilvert, D. Is chorismate Mutase a Prototypic Entropy Trap? - Activation Parameters for the Bacillus Subtilis Enzyme. Tetrahedron Lett. 1996, 37, 2691–2694.

29



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(41) Ishida, T.; Fedorov, D. G.; Kitaura, K. All Electron Quantum Chemical Calculation of the Entire Enzyme System Confirms a Collective Catalytic Device in the Chorismate Mutase Reaction. J. Phys. Chem. B. 2006, 110, 1457–1463. (42) Kienh¨ofer, A.; Kast, P.; Hilvert, D. Selective Stabilization of the Chorismate Mutase Transition State by a Positively Charged Hydrogen Bond Donor. J. Am. Chem. Soc. 2003, 125, 3206–3207. (43) Lyne, P. D.; Mulholland, A. J.; Richards, W. G. Insights into Chorismate Mutase Catalysis from a Combined QM/MM Simulation of the Enzyme Reaction. J. Am. Chem. Soc. 1995, 117, 11345–11350. (44) Claeyssens, F.; Harvey, J. N.; Manby, F. R.; Mata, R. A.; Mulholland, A. J.; Ranaghan, K. E.; Sch¨ utz, M.; Thiel, S.; Thiel, W.; Werner, H.-J. High-Accuracy Computation of Reaction Barriers in Enzymes. Angew. Chem. 2006, 118, 7010–7013. (45) Claeyssens, F.; Ranaghan, K. E.; Manby, F. R.; Harvey, J. N.; Mulholland, A. J. Multiple High-Level QM/MM Reaction Paths Demonstrate Transition-State Stabilization in Chorismate Mutase: Correlation of Barrier Height with Transition-State Stabilization. Chem. Comm. 2005, 5068. (46) Fedorov, D. G.; Kitaura, K. Pair Interaction Energy Decomposition Analysis. J. Comput. Chem. 2006, 28, 222–237. (47) Morokuma, K. Molecular Orbital Studies of Hydrogen Bonds. III. C=O···H-O Hydrogen Bond in H2CO···H2O and H2CO···2H2O. J. Chem. Phys. 1971, 55, 1236. (48) Ziegler, T.; Rauk, A. On the Calculation of Bonding Energies by the Hartree Fock Slater Method. Theoret. Chim. Acta 1977, 46, 1–10. (49) Phipps, M. J. S.; Fox, T.; Tautermann, C. S.; Skylaris, C.-K. Energy Decomposition

30


Page 30 of 42

Page 31 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Analysis Approaches and Their Evaluation on Prototypical Protein–Drug Interaction Patterns. Chem. Soc. Rev. 2015, 44, 3177–3211. (50) Heifetz, A.; Chudyk, E. I.; Gleave, L.; Aldeghi, M.; Cherezov, V.; Fedorov, D. G.; Biggin, P. C.; Bodkin, M. J. The Fragment Molecular Orbital Method Reveals New Insight into the Chemical Nature of GPCR–Ligand Interactions. J. Chem. Inf. Model 2016, 56, 159–172. (51) Nakano, T.; Kaminuma, T.; Sato, T.; Fukuzawa, K.; Akiyama, Y.; Uebayasi, M.; Kitaura, K. Fragment Molecular Orbital Method: Use of Approximate Electrostatic Potential. Chem. Phys. Lett. 2002, 351, 475–480. (52) Fedorov, D. G.; Kitaura, K. The Importance of Three-Body Terms in the Fragment Molecular Orbital Method. J. Chem. Phys. 2004, 120, 6832. (53) Fedorov, D. G.; Kitaura, K. The Three-Body Fragment Molecular Orbital Method for Accurate Calculations of Large Systems. Chem. Phys. Lett. 2006, 433, 182–187. (54) Fedorov, D. G.; Slipchenko, L. V.; Kitaura, K. Systematic Study of the Embedding Potential Description in the Fragment Molecular Orbital Method †. J. Phys. Chem. A 2010, 114, 8742–8753. (55) Fedorov, D. G.; Jensen, J. H.; Deka, R. C.; Kitaura, K. Covalent Bond Fragmentation Suitable To Describe Solids in the Fragment Molecular Orbital Method. J. Phys. Chem. A 2008, 112, 11808–11816. (56) Fedorov, D.; Kitaura, K. Second Order Møller-Plesset Perturbation Theory Based Upon the Fragment Molecular Orbital Method. J. Chem. Phys. 2004, 121, 2483. (57) Fedorov, D.; Ishimura, K.; Ishida, T.; Kitaura, K.; Pulay, P.; Nagase, S. Accuracy of the Three-Body Fragment Molecular Orbital Method Applied to Møller–Plesset Perturbation Theory. J. Comput. Chem. 2007, 28, 1476–1484. 31



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(58) Fletcher, G. D.; Fedorov, D. G.; Pruitt, S. R.; Windus, T. L.; Gordon, M. S. LargeScale MP2 Calculations on the Blue Gene Architecture Using the Fragment Molecular Orbital Method. J. Chem. Theory Comput. 2012, 8, 75–79. (59) Pruitt, S. R.; Nakata, H.; Nagata, T.; Mayes, M.; Alexeev, Y.; Fletcher, G.; Fedorov, D. G.; Kitaura, K.; Gordon, M. S. Importance of Three-Body Interactions in Molecular Dynamics Simulations of Water Demonstrated with the Fragment Molecular Orbital Method. J. Chem. Theory Comput. 2016, 12, 1423–1435. (60) Fedorov, D. G.; Ishida, T.; Kitaura, K. Multilayer Formulation of the Fragment Molecular Orbital Method (FMO). J. Phys. Chem. A 2005, 109, 2638–2646. (61) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S. et al. General Atomic and Molecular Electronic Structure System. J. Comput. Chem. 1993, 14, 1347–1363. (62) Hehre, W. J. Self-Consistent Molecular-Orbital Methods. I. Use of Gaussian Expansions of Slater-Type Atomic Orbitals. J. Chem. Phys. 1969, 51, 2657. (63) Hehre, W. J.; Ditchfield, R.; Pople, J. A. Self Consistent Molecular Orbital Methods. XII. Further Extensions of Gaussian–Type Basis Sets for Use in Molecular Orbital Studies of Organic Molecules. J. Chem. Phys. 1972, 56, 2257– 2261. (64) Dunning, T. H. Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron Through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007. (65) Woon, D. E.; Dunning, T. H. Gaussian Basis Sets for Use in Correlated Molecular Calculations. III. The Atoms Aluminum Through Argon. J. Chem. Phys. 1993, 98, 1358. (66) Jensen, F. Unifying General and Segmented Contracted Basis Sets. Segmented Polarization Consistent Basis Sets. J. Chem. Theory Comput. 2014, 10, 1074–1085. 32


Page 32 of 42

Page 33 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(67) Jensen, F. Segmented Contracted Basis Sets Optimized for Nuclear Magnetic Shielding. J. Chem. Theory Comput. 2015, 11, 132–138. (68) Steinmann, C.; Ibsen, M. W.; Hansen, A. S.; Jensen, J. H. FragIt: A Tool to Prepare Input Files for Fragment Based Quantum Chemical Calculations. PLoS ONE 2012, 7, e44480. (69) Edmiston, C.; Ruedenberg, K. Localized Atomic and Molecular Orbitals. Rev. Mod. Phys. 1963, 35, 457–464. (70) Fedorov, D. G.; Kitaura, K. Use of an Auxiliary Basis Set to Describe the Polarization in the Fragment Molecular Orbital Method. Chem. Phys. Lett. 2014, 597, 99 – 105. (71) Ranaghan, K. E.; Ridder, L.; Szefczyk, B.; Sokalski, W. A.; Hermann, J. C.; Mulholland, A. J. Insights Into Enzyme Catalysis from QM/MM Modelling: Transition State Stabilization in Chorismate Mutase. Mol. Phys. 2003, 101, 2695–2714. (72) Ranaghan, K. E.; Mulholland, A. J. Simulating Enzyme Reactivity; Royal Society of Chemistry (RSC), pp 375–403.

33




Page 34 of 42

Page 35 of 42



A

10 5 0 −5 −10 −15

20

6-31G(d) cc-pVDZ

pcseg-1

−20 −2.0 −1.5 −1.0 −0.5 0.0 0.5 1.0 1.5 ˚ Reaction Coordinate [A]

15

Page 36 of 42

FMO2-MP2/AFO

3

B

10 5 0 −5 −10 −15 −20 −2.0 −1.5 −1.0 −0.5 0.0 0.5 1.0 1.5 ACS Paragon Plus Environment ˚ Reaction Coordinate [A]

Barrier Difference [kcal/mol]

15

FMO2-MP2/HOP Reaction Barrier [kcal/mol]

Reaction Barrier [kcal/mol]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

20


C

2 1 0 −1 −2 −3 −2.0 −1.5 −1.0 −0.5 0.0 0.5 1.0 1.5 ˚ Reaction Coordinate [A]

Error From FMO2-MP2/6-31G(d)

A

2 1 0 −1 −2

3

6-31G(d):6-31G(d) STO-3G:6-31G(d) pcseg-0:6-31G(d)


Error From FMO2-MP2/cc-pVDZ

B

2 1 0 −1 −2

3

6-31G(d):cc-pVDZ STO-3G:cc-pVDZ cc-pVDZ:cc-pVDZ

−3 −2.0 −1.5 −1.0 −0.5 0.0 0.5 1.0 1.5 ACS Paragon Plus Environment ˚ Reaction Coordinate [A]

Barrier Error [kcal/mol]


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

3



Page 37 of 42

Error From FMO2-MP2/pcseg-1

C

2 1 0 −1 −2

pcseg-0:pcseg-1 pcseg-1:pcseg-1




Page 38 of 42

Page 39 of 42

E es

20

Arg7

B 1.5

10

Interaction Energy [kcal/mol]


15

E disp

2.0 Arg90 Arg63

A

5 0 −5 −10 −15

1.0 0.5 0.0 −0.5 −1.0 −1.5

−20 −2.0

−1.5

−1.0 −0.5 0.0 0.5 1.0 ˚ Reaction Coordinate [A]

1.5

−2.0 −2.0

E ct+mix

2.0

−1.5


1.5

E ex

10

C

D

1.5 1.0



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.5 0.0 −0.5 −1.0

5

0

−5

−1.5 −2.0 −2.0

−1.5


1.5

−10 −2.0


−1.5


1.5



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

Page 40 of 42

25

E es E ex E disp E ct

20 15 10 5 0 −5 −10

1

2

3 ACS Paragon 4 Plus 5 Environment 6 7

8

Avg

Page 41 of 42

1 2 3 4 5 6 7 8 9 10 11 12 13 O14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48


O-

O-

O O

4

O

3

3 2O

1 O O-

O

2

O

1

2 OH

4 4

O

3

-O

O-

O-

1

HO OH

1

2


3



Page 42 of 42

Mapping Interaction Energies in Chorismate Mutase with the Fragment

Recommend Documents