Catalytic Effect of Aqueous Solution in Water-Assisted Proton-Transfer

Mar 8, 2018 - College of Physics and Electronics, Shandong Normal University, Jinan ... The transition-state dipole moment is the largest along the re...
0 downloads 0 Views 6MB Size
Article Cite This: J. Phys. Chem. B 2018, 122, 3124−3132

pubs.acs.org/JPCB

Catalytic Effect of Aqueous Solution in Water-Assisted ProtonTransfer Mechanism of 8‑Hydroxy Guanine Radical Peng Liu, Chen Li, Shengyu Wang, and Dunyou Wang* College of Physics and Electronics, Shandong Normal University, Jinan 250014, China

Downloaded via KAOHSIUNG MEDICAL UNIV on July 9, 2018 at 07:26:59 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Water-assisted proton-transfer process is a key step in guanine damage reaction by hydroxyl radical in aqueous solution. In this article, we quantitatively determine the solvent effect in waterassisted proton-transfer mechanism of 8-hydroxy guanine radical using combined quantum mechanics and molecular mechanism with an explicit solvation model. Atomic-level reaction pathway was mapped, which shows a synchronized two-proton-transfer mechanism between the assistant water molecule and 8-hydroxy guanine radical. The transition-state dipole moment is the largest along the reaction pathway, which electrostatically stabilizes the proton-transfer transition-state complex. The free-energy reaction barrier for this waterassisted proton-transfer reaction was calculated at 19.2 kcal/mol with the density functional theory/M08-SO/cc-pVTZ+/molecular mechanics level of theory. The solvent effect not only has a big impact on geometries, but also dramatically changes the energetics along the reaction pathway. Among the solvent effect contributions to the transition state, the solvent energy contribution is −28.5 kcal/mol and the polarization effect contribution is 19.9 kcal/mol. In total, the solvent effect contributes −8.6 kcal/mol to the free-energy barrier height, which means that the presence of aqueous solution has a catalytic effect on the reaction mechanism and enhances the proton-transfer reactivity in aqueous solution.

I. INTRODUCTION

attacks the substrate CCl4. Furthermore, a solvent can serve as a solvent as well as a catalyst at the same time.11−13 Water-assisted proton-transfer mechanisms have been observed in the reactions with RNA and DNA bases, including uracil,14 cytosine,15−18 thymine,18,19 adenine,18,20 and guanine.18,21−25 Quantum chemical study of intramolecular proton transfer for monohydrated guanine complexes found that the reaction barrier height is approximately two times lower for the tautomeric oxo-hydro reaction compared to non-water-assisted processes.21 Ab initio study at the density functional theory (DFT) and MP2 levels of theory of the intramolecular proton transfer in C8-oxidative guanine22 showed that the assistance of a water molecule in the proton-transfer process greatly reduces the reaction barrier; furthermore, the presence of oxygen at C8 position lowers the activation energy by about 1 kcal/mol. First-principles study using DFT and a continuum solvent model (integral equation formalism-polarizable continuum model (IEF-PCM))23 found that a water-molecule-assisted proton-transfer process of tautomerization in 8-oxo-7,8dehydropurine can significantly lower the activation Gibbs energy to less than 10 kcal/mol in aqueous solution. As for 8hydroxy guanine radical (8-OHGrad), ab initio multiconfigurational MCSCF/DZP/MRPT2 study25 showed the activation

The presence of a solvent can change the energetics of the chemical reactions dramatically from their corresponding ones in gas phase. Combined multilevel quantum mechanics (QM) and molecular dynamics studies1−4 show that the presence of an aqueous solution weakens the interactions between a nucleophile and a substrate at its transition state (TS) for the bimolecular nucleophilic substitution (SN2) reaction, which substantially raises the reaction barrier height compared to the one in gas phase. For example, for the normal Walden-inversion mechanism of the F− + CH3Cl reaction mechanism, the reaction barrier in gas phase is ∼ 6.9 kcal/mol;5 however, it is substantially increased to 26.9 kcal/mol6 because of the presence of the aqueous solution. As a result, the reaction rate constant in gas phase is about 20 orders of magnitude larger than that in solution phase for this SN2 reaction. A solvent can also change the gas-phase reaction mechanism by participating in a chemical reaction directly. For example, for the tautomerization of formamide to formamidic acid reaction,7−9 it is found that the water molecule assists tautomerization of formamide in aqueous solution and lowers the barrier height compared to the one in gas phase. Another study of the water-assisted SN2 reaction of OH− + CCl4 in aqueous solution10 found that a water molecule assists to form a new OH− in the favorable back-side attack conformation via the proton-transfer process and then the newly formed OH− © 2018 American Chemical Society

Received: October 9, 2017 Revised: March 3, 2018 Published: March 8, 2018 3124

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

where VQM has the same expression as in the gas phase, representing the potential of QM region, and VMM describes the classical molecular mechanical energy of the MM region. The middle term, VQM/MM, is the interaction energy between the QM and MM subsystems, including the bonded interactions, electrostatic interactions, and van der Waals interactions. Moreover, the electrostatic interactions between the solute electron density ρ and the solvent classical charge ZI can be approximated as32

barrier of water-assisted proton transfer in gas phase is 22.4 kcal/mol, which is about 26.3 kcal/mol less than the direct reaction process without the assistance of one water molecule. In aqueous solution, a DFT/B3LYP theory with the IEFpolarizable continuum solvation model24 found that the waterassisted proton transfer of 8-OHGrad with one water molecule is at 18.6 kcal/mol. Although the authors characterized the stationary points along the reaction pathway of the waterassisted proton-transfer process, the detailed, atomic-level evolution of the reaction pathway was not available, and the solvent effects and the quantitative contributions of the solvent to the reaction barrier are still unknown. It is well known that the proton-transfer process has to be assisted by water molecules’ dipole moment distribution from reorganization of the first and second solvation shells; thus, to correctly reproduce the microscopic rearrangements of water molecules, one should treat solution water molecules explicitly. However, to date, there have been no studies on water-assisted proton transfer of 8-OHGrad in aqueous solution to analyze the quantitative contribution of the solvent effect to the energetics along the reaction pathway; therefore, the catalytic effect of the aqueous solution on the activation barrier was unknown. In addition, there has been no detailed reaction pathway mapped. Thus, in this paper, using an explicit SPC/E solvent model,26 we combined quantum mechanics and molecular mechanics (MM) theory1,27,28 to study the waterassisted intramolecular proton transfer of 8-OHGrad in water. Our purpose is (i) not only to investigate the solvent effects on this reaction mechanism in aqueous solution, (ii) but also to calculate the quantitative contribution of the solvent effect to the potential of mean force (PMF), (iii) to map a detailed, atomic-level, water-assisted proton-transfer mechanism in aqueous solution, and (iv) to calculate the potential of mean force (PMF) and free-energy activation barrier of the title reaction in aqueous solution.

Velectrostatic =

I

=

∑ i,I

ZI Q i |R I − ri|

ZIρ(r′) dr′ |R I − r′| = VESP(r, R, Q ) (2)

Here, ESP charges Qi stands for the QM electronic density, which was fitted from the ab initio calculation of the electrostatic potential surface of the QM region. The potential of mean force under the DFT/MM level of theory can be obtained as27 DFT DFT ← ESP DFT ← ESP ESP ΔWAB = (ΔWAA + ΔWBB ) + ΔWAB

(3)

where the last term is the solvent contribution to the PMF, calculated by statistical samplings with classical ESP/MM description. The first term in parentheses denotes the internal energy of the solute region under DFT level of theory. First, we embedded the estimated reaction complex (RC) in water, optimizing the QM and MM subsystems with a multiregional optimization protocol. Then, we fixed the solute region and equilibrated the solvent MM region using molecular dynamics simulation for 120 ps at 298 K, where the fixed QM region was represented by the ESP charges obtained from the previous optimization step. Second, the whole system was optimized again and the optimized structure was used to search the product complex (PC) with a bond-breaking and bondformation procedure.2 This PC was optimized and equilibrated as we did to the initial RC. Third, the transition state (TS) was identified based on the optimized RC and PC, and was confirmed using a numerical frequency calculation with one imaginary frequency. Subsequently, we obtained the final RC and PC by optimizing displacements of the TS along the imaginary frequency mode. Finally, using the final obtained RC and PC, the reaction pathway was mapped using the nudged elastic band (NEB) approach.33 Molecular dynamics simulation was used to equilibrate the solvent for 120 ps along the whole reaction pathway, and the whole NEB reaction pathway was then optimized. The last step was repeated until the final NEB reaction pathway was converged. Finally, the PMF was calculated with the DFT/MM level of theory according to eq 3.

II. METHODS Using combined quantum mechanics and molecular mechanics theory,1,27,28 we investigated the water-assisted intramolecular proton transfer of 8-OHGrad in water, where the 8-OHGrad complex with one water molecule was treated as the QM region with DFT and effective electrostatic potential (ESP) levels of theory during different stages of calculation, and the aqueous solution was treated as the classical MM region. The 8OHGrad·H2O complex was embedded in a 39.9 Å cubic box consisting of 2137 water molecules described with an explicit SPC/E water model.26 As Truhlar and co-workers pointed out, the DFT/M08-SO/cc-pVTZ+ calculation gives comparable accuracy in terms of reaction barrier height as the CCSD(T)(full)/aug-cc-pCV(T + d)Z calculation.29 Therefore, in this paper, the DFT/M08-SO/cc-pVTZ+ level of theory was employed to treat the QM region to calculate PMFs along the reaction pathway. In addition, the solute interacts with water molecules via bonded interactions, electrostatic interactions, and van der Waals interactions, and the van der Waals parameters of the QM region were obtained from standard Amber force field.30 The calculations of PMFs were performed using the NWChem computational chemistry package.31 The potential energy of the whole reaction system can be written as Vpotential = VQM + VQM/MM + VMM

∑∫

III. RESULTS AND DISCUSSION III.I. Stationary Points along the Reaction Pathway. First of all, the •OH attacks C8 site of guanine and results in the reactant (8-OHGrad) for the proton-transfer process. The •OH is attached to the C8 site in guanine, as shown in Figure 1. In aqueous solution, the hydroxyl radical, O10H10, is attached to C8 with a bond distance of 1.402 Å, and the bending angle of ∠C8O10H10 is 108.1°; the dihedral angle of ∠C5N7C8O10 is 123.3°; and the dihedral angle with respect the proton H10 ∠N7C8O10H10 is 15.5°. The location of O10H10 radical in gas phase is similar to the one in the solution phase. The distance of C8O10 is 1.382 Å, and the bending angle of

(1) 3125

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

between the hydroxyl located at the C8 site and its adjacent water molecules, with an average distance of 1.939 Å. There are also two weak hydrogen bonds formed with N7, with an average distance of 2.653 Å. Thus, the interactions between the solute and solvent, especially the hydrogen bonds, affect the reactant geometry in aqueous solution. And solvent effect here has a much bigger impact on the torsion angles than on the bond distances and bending angles. Figure 2 shows the structures of transition state in gas phase25 and in aqueous solution, where the latter was

Figure 1. Solvent effects on geometry. Comparison of reactant complexes for water-assisted proton-transfer process of 8-OHGrad in (a) gas phase and (b) aqueous solution. The indicated distances are in angstroms (bond lengths in black and hydrogen bonds in blue).

∠C8O10H10 is 108.8°; the dihedral angle of ∠C5N7C8O10 is 110.9°; and the dihedral angle with respect the proton H10 ∠N7C8O10H10 is 58.0°. (The structures of the stationary points are provided in Table S1 of the Supporting Information.) Comparisons of stationary points for water-assisted intramolecular proton-transfer process of 8-OHGrad in gas phase and in aqueous solution are shown in Figure 1−3, where the stationary points in gas phase are from Chaban et al.25 The optimized structures of reactant complex in gas phase and in solution phase are compared in Figure 1. The distance of H10O11 is 1.823 Å in aqueous solution, which is about 0.236 Å shorter than that in gas phase; the distance of N7H11 is 2.045 Å in aqueous solution, which is about 0.162 Å shorter than that in gas phase. The biggest differences exist among dihedral angles between in gas phase and in solution phase: the dihedral angles of ∠N7C8O10H10, ∠C8O10H10O11, ∠O10H10O11H11, ∠H10O11H11N7, ∠O11H11N7C8, and ∠H11N7C8O10 are 58.0, −26.2, −10.9, 11.8, 13.4, and −51.4°, respectively, in gas phase, whereas they are 15.5, 5.3, −2.4, −17.7, 35.4, and −25.3°, respectively, in aqueous solution, with a mean unsigned difference of 26.7° with respect to those in gas phase. We found that there are hydrogen bonds formed between N7, C8, O10, O11, H11 sites and the adjacent water molecules: there are two strong hydrogen bonds formed with the H11O11H12 water molecule with an average distance of 1.746 Å; one strong and one weak hydrogen bond formed

Figure 2. Comparison of transition states for water-assisted protontransfer process of 8-OHGrad in (a) gas phase and (b) aqueous solution. The indicated distances are in angstroms (bond lengths in black and hydrogen bonds in blue).

confirmed by numerical frequency calculations with a single imaginary frequency of 793.7i cm−1. In transition state, the water molecule and the hydroxyl radical, O10H10, form a closed ring with N7C8 as H11 in the water molecule attached to N7 and O11 in the hydroxyl radical attached to H10. The biggest difference between these two transition states still remains in their dihedral angles. In gas phase, the dihedral angles of ∠N7C8O10H10, ∠C8O10H10O11, ∠O10H10O11H11, ∠H10O11H11N7, ∠O11H11N7C8, and ∠H11N7C8O10 are 32.5, −9.3, −12.7, 8.2, 11.9, and −34.4°, respectively, whereas in solution phase, they are 29.1, −26.0, 10.9, −11.2, 28.1, and −33.6° respectively, with a mean unsigned difference of 13.4° with respect to those in gas phase. Again, we note that four strong hydrogen bonds are formed, one formed with O11 and three formed with O10, with an average distance of 1.673 Å; meanwhile, three weak hydrogen bonds, one with C8H8 and two with N7, occur with an average 3126

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

bond distances, which is different from our previous investigation of SN2 reactions1−4 that the solvent effect can have a big influence on not only the angles but also the bond lengths in aqueous solution. However, the current reaction system has a compact ring structure and the solvent effect is not strong enough to affect the strong bond interactions in the stationary points, whereas the more weaker torsion interactions can be easily affected by the solvent effect, especially by the hydrogen bonds formed in their hydrogen shells. Therefore, for the current system, the biggest difference of the conformations of the stationary points lies in their corresponding torsion angles between gas phase and solution phase. Mulliken spin density distributions of the reactant complex, transition state, and product state were also calculated using the Gaussian program.34 We note that the hybrid density functional M08-SO is not included in the Gaussian program, instead we chose the IEF-PCM/M06-2X/aug-cc-pVTZ method to calculate the spin densities, as Truhlar and co-workers pointed out that the M06-2X/aug-cc-pVTZ level of theory produces very close results to the M08-SO/cc-pVTZ+.29 Figure 4 shows the spin density distributions of the reactant complex, transition state, and product complex. As the OH radical is attached to the C8 site to form the reactant complex, 8-hydroxy guanine radical (8-OHGrad), the localized spin density of a separate

distance of 2.711 Å. Therefore, the big changes of torsion angles here compared to gas phase are strongly influenced by the hydrogen bonds formed between the solute and solvent. The product complexes, both in the gas25 and solution phases, are plotted in Figure 3. In the product state, H11,

Figure 3. Comparison of product complexes for water-assisted protontransfer process of 8-OHGrad in (a) gas phase and (b) aqueous solution. The indicated distances are in angstroms (bond lengths in black and hydrogen bonds in blue).

originated from the assistant water molecule, has been transferred to the N7 site, whereas H10, from hydroxyl on the C8 site, has been transferred to the O11 and formed a new water molecule, H10O11H12. Thus, the proton-transfer process has finished with the assistance of one water molecule in aqueous solution. At this stage, the proton H10 has been broken up from the hydroxyl radical and been transferred to O11. However, the biggest difference still exists in their dihedral angles: ∠N7C8O10H10, ∠C8O10H10O11, ∠O10H10O11H11, ∠H10O11H11N7, ∠O11H11N7C8, and ∠H11N7C8O10 are 35.0, −15.5, −4.5, 2.9, 15.8, and −38.5°, respectively, in gas phase, whereas they are 41.5, −10.0, −13.0, 3.8, 25.2, and −53.5°, respectively, in aqueous solution, with a mean unsigned difference of 7.6° between them. Two strong hydrogen bonds formed around the newly formed water molecule with an average distance of 1.812 Å; there are also two strong hydrogen bonds formed at O10 with an average distance of 1.648 Å, and five weak hydrogen bonds, three formed at N7H11 and two at C8H8, with an average distance of 2.652 Å. As a result, the solvent effect has a tremendous impact on the torsion angles, much more than on the bending angles and the

Figure 4. Spin density distribution of the (a) reactant complex, (b) transition state, and (c) product complex of QM solute obtained at the IEF-PCM/M06-2X/aug-cc-pVTZ level of theory. 3127

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

complex, which is the process of the breaking of old OH bond in the assistant water molecule and the forming of H11N7 bond. This is a concerted two-proton-transfer process. Furthermore, there are also big changes in their bending and dihedral angles during the proton-transfer process. For instance, the angle of ∠C8O10H10 decreases from 108.1° at reactant complex to 105.9° at transition state, then to 101.8° at product complex; and the angle of ∠O10H10O11 decreases from 164.9° at reactant complex to 163.8° at transition state, then to 151.3° at product complex. In addition, the dihedral angle of ∠N7C8O10H10 increases from 15.5° at reactant complex to 29.1° at transition state, then to 41.5° at product complex; the dihedral angle of ∠H10O11H11N7 increases from −17.7° at reactant complex to −11.2° at transition state, then to 3.8° at product complex; the dihedral angle of ∠O11H11N7C8 decreases from 35.4° at reactant complex to 28.1° at transition state, then to 25.2° at product complex; and the dihedral angle of ∠H11N7C8O10 decreases from −25.3° at reactant complex to −33.6° at transition state, then to −53.5° at product complex. In short, the proton-transfer process is assisted with one water molecule, which is a synchronized two-proton-transfer process: one proton is transferred from the hydroxyl to the assistant water molecule and the other from the assistant water molecule to the N7 site. This whole proton-transfer process involves two bonds breaking and forming concertedly. The dipole moment evolution of the solute subsystem along the NEB reaction pathway is plotted in Figure 6. It shows that at the transition state of point (6) the 8-OHGrad with the assistant water molecule has the largest dipole moment. This means that at the transition state the solute subsystem has the strongest electrostatic interaction with the surrounding water molecules in its solvation shells. The strongest electrostatic interaction here stabilizes the transition state, which indicates that the activation energy will be lowered due to the strongest interaction between the solute and solvent at the transition state. In addition, the dipole moment of the reactant complex is smaller than the product complex because the reactant complex is more compact than the product complex. III.III. Potentials of Mean Force and Solution Contributions. We plotted the potential of mean force (PMF) as well as the solvent contribution along the NEB reaction pathway with the reactant as a reference point in Figure 7. It shows that the free-energy reaction barrier is 19.2 kcal/mol and the reaction energy is 5.5 kcal/mol. Moreover, the water solvent contributes −28.5 kcal/mol to the TS and −27.3 kcal/ mol to the product complex, which indicates that the presence of solvent significantly reduces the reaction barrier height. This is consistent with our prediction on the analysis of the dipole moment that the largest electrostatic interaction with the surrounding water molecules at the transition state lowers the barrier height. The perturbation of the solute wave function due to the presence of the solvent leads to the contribution of polarization effect, which is shown in Figure 8 by comparing the gas-phase reaction pathway and the solute internal energy along the NEB reaction pathway. Here, the gas-phase energy is obtained with the same 10 snapshots along the PMF reaction pathway without the solvent contribution and solute−solvent interactions; the solute internal energy is obtained excluding the solvent energy contribution. The comparison shows that the polarization effect contributes 32.2 kcal/mol to reactant state, 52.1 kcal/mol to transition state, and 52.9 kcal/mol to product state. Thus, the net polarization effect contributes 19.9

O10H10 radical has been delocalized and widely distributed in the 8-OHGrad radical as O10 only has a small spin density of 0.040, whereas N7 has the largest spin density of 0.507 (Table 1). On the one hand, as the reaction proceeds from reactant Table 1. Mulliken Spin density of the Reactant Complex (RC), Transition State (TS), and Product Complex (PC) of QM Solute Obtained at the IEF-PCM/M06-2X/aug-ccpVTZ Level of Theory H9 N9 C8 H8 N7 C5 C6 O6 N1 H1 C2 N2 1H2 2H2 N3 C4 O10 H10 H11 O11 H12

RC

TS

PC

−0.00387 0.05800 0.00011 0.01726 0.50741 0.01026 0.08354 0.11623 −0.01468 −0.00173 0.04924 0.07398 −0.00627 −0.00620 −0.03058 0.10939 0.03997 −0.00683 −0.00156 0.00650 −0.00016

−0.00688 0.08265 −0.03917 0.01472 0.49693 0.24035 0.03501 0.13259 −0.00190 −0.00202 0.02484 0.06377 −0.00434 −0.00404 −0.00830 −0.00169 0.02838 −0.00622 −0.04529 −0.00003 0.00063

−0.01228 0.09599 −0.00507 0.00624 0.47369 0.28262 0.02806 0.13575 −0.00144 −0.00150 0.03264 0.05454 −0.00329 −0.00210 −0.01176 −0.00879 0.01723 −0.00457 −0.07706 0.00126 −0.00015

complex to transition state, then to product complex, the spin density of O10 has been further decreased from 0.040 to 0.028, then to 0.017, and H10 from −0.007 to −0.006, then to −0.005. On the other hand, because H11 is transferred to the N7 site, it gains spin density from −0.002 at reactant complex to −0.045 at transition state, then to −0.078 at product state, whereas the spin density of N7 decreases from 0.507 to 0.497, then to 0.474. III.II. Detailed, Atomic-Level Reaction Mechanism. Ten snapshots along the NEB33 reaction pathway are displayed in Figure 5 to examine the atomic-level evolution of the reaction in detail. Snapshot (1) shows the reactant structure, snapshot (6) shows the transition state, and snapshot (10) shows the product structure. First, from snapshot (1) to snapshot (10), the charges distributed on H10 are 0.519, 0.449, 0.498, 0.512, 0.549, 0.559, 0.530, 0.544, 0.621, and 0.481; the ones distributed on H11 are 0.381, 0.319, 0.408, 0.361, 0.470, 0.504, 0.442, 0.489, 0.373, and 0.48. Thus, the reaction mechanism is characterized by a synchronized two-protontransfer process: one proton transfer involves H10 transferred to the assistant water molecule, whereas the other involves the transfer of H11 in the assistant water molecule to N7 site concertedly. As can been seen, the H10O11 bond distance is 1.823 Å at reactant complex, which decreases to 1.164 Å at transition state and further decreases to 0.980 Å at product complex, which is the process of breaking of the O10H10 bond in the 8-OHGrad and forming of the new OH bond in the assistant water molecule. At the meantime, the H11N7 bond distance is 2.045 Å at reactant complex, decreases to 1.132 Å at transition state, and further decreases to 1.018 Å at product 3128

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

Figure 5. Structures of 10 snapshots along the NEB reaction pathway for the water-assisted proton transfer of 8-OHGrad in aqueous solution. Snapshot (1) is the structure of reactant complex; snapshot (6), transition state; and snapshot (10), product complex. The indicated distances are in angstroms. 3129

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

presence of aqueous solution enhances the reactivity. In other words, the aqueous solution has a catalytic effect on this reaction mechanism by reducing the reaction barrier about 8.6 kcal/mol. This is also confirmed by comparison with the reaction barrier height in gas phase. In gas phase, the title reaction has a barrier height at 22.4 kcal/mol,25 which is lowered to 19.2 kcal/mol here in aqueous solution. Therefore, the presence of aqueous solution indeed enhances the protontransfer reactivity. Our calculated barrier height at 19.2 kcal/mol has a very good agreement with the one at 18.6 kcal/mol by Munk et al.24 using a continuum solvation model with IEF-PCM/B3LYP/ aug-cc-pVTZ//B3LYP/6-31G(d) level of theory. Usually, a system in solution treated with implicit and explicit solvent models would produce quite different results. For example, as Yang’s group35 compared the reactant structure of the Cl− + CH3Cl in solution optimized with implicit and explicit solvent models, they found that the reactant structure optimized using implicit solvent model leads to a wrong symmetry conformation. Therefore, they pointed out that there are “risks of using conformations optimized in the gas phase or with a continuum solvent model for the study of solution or enzyme reactions”. In addition, the ring-opening process of guanine damage by hydroxyl radical36 with an explicit SPC/E water model shows that the reaction barrier height is 31.6 kcal/mol at the CCSD(T)/MM level of theory, which is much higher than the one at 19.5 kcal/mol obtained with the IEF-PCM implicit water model at DFT/B3LYP level of theory.24 Thus, the good agreement of the barrier height in this work between different solvent models might be a coincidence caused by the computational error produced from the implicit solvent model and the quantum level of theory used there. Our calculated free reaction energy, 5.5 kcal/mol, is about 2.5 kcal/ mol lower than the one at 8.0 kcal/mol calculated by Munk et al. using a continuum solvation model. Because there are no experimental data of this reaction to compare with, here we used the method (the implicit, ICEPCM solvation model in Gaussian program34) employed by Munk et al.24 to try to reproduce our explicit water model results to show the consistency between the two methods. On the basis of the same QM structures of the stationary points in aqueous solution from our explicit SPC/E water model results, the IEF-PCM/M06-2X/aug-cc-pVTZ calculation gives the barrier height at 19.3 kcal/mol and the free reaction energy at 7.1 kcal/mol, which are consistent with our explicit water model results at 19.5 and 5.5 kcal/mol, respectively, at the DFT/M08-SO/cc-pVTZ+ level of theory.

Figure 6. Evolution of the dipole moment of the solute subsystem along the NEB reaction pathway for the water-assisted proton transfer of 8-OHGrad in aqueous solution.

Figure 7. Potential of mean force and solvent contributions along the NEB reaction pathway. The potential of mean force was calculated at the DFT/MM level of theory with the reactant state as the energy reference point.

IV. SUMMARY AND CONCLUSIONS Combining quantum mechanics with molecular mechanics theories, we investigated the water-assisted intramolecular proton transfer of 8-OHGrad in aqueous solution. An explicit SPC/E solvation model26 was used to describe the aqueous solution to describe the local rearrangements of water molecules. The stationary points in comparison to the corresponding ones in gas phase show that the aqueous solution affects the dihedral angles the most. Detailed, atomiclevel reaction mechanism along the reaction pathway shows that the reaction mechanism is a synchronized two-protontransfer process: on the one hand, the assistant water molecule accepts proton H10 from the hydroxyl attached at C8 site; on the other hand, it transfers its proton H11 to the N7 site of 8OHGrad.

Figure 8. Comparison between gas phase and solute internal energies along the NEB reaction pathway under the DFT/MM level of theory using the gas-phase energy of the reactant complex as a reference point.

kcal/mol to transition state and 20.7 kcal/mol to product state. In total, the aqueous solution contributes −8.6 kcal/mol to free-energy barrier height and −6.6 kcal/mol to reaction energy. So, the solvent energy and polarization effect have opposite contributions to the PMF: the former has an effect on reducing the barrier height, and the latter raises the barrier. Nonetheless, the total solvent effect has contributed −8.6 kcal/ mol to the free-energy barrier height, which means that the 3130

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B

(8) Markova, N.; Enchev, V. Water-assisted Proton Transfer in Formamide, Thioformamide and Selenoformamide. J. Mol. Struct.: THEOCHEM 2004, 679, 195−205. (9) Liang, W.; Li, H.; Hu, X.; Han, S. Proton Transfer of Formamide + nH2O (n = 0−3): Protective and Assistant Effect of the Water Molecule. J. Phys. Chem. A 2004, 108, 10219−10224. (10) Chen, J.; Yin, H.; Wang, D.; Valiev, M. Water Assisted Reaction Mechanism of OH− with CCl4 in Aqueous Solution − Hybrid Quantum Mechanical and Molecular Mechanics Investigation. Chem. Phys. Lett. 2013, 559, 30−34. (11) Cole, A. C.; Jensen, J. L.; Ntai, I.; Tran, K. L. T.; Weaver, K. J.; Forbes, D. C.; Davis, J. H. Novel Brønsted Acidic Ionic Liquids and Their Use as Dual Solvent− Catalysts. J. Am. Chem. Soc. 2002, 124, 5962−5963. (12) Caló, V.; Nacci, A.; Monopoli, A.; Fanizzi, A. Cyclic Carbonate Formation from Carbon Dioxide and Oxiranes in Tetrabutylammonium Halides as Solvents and Catalysts. Org. Lett. 2002, 4, 2561−2563. (13) Kim, D. W.; Ahn, D.-S.; Oh, Y.-H.; Lee, S.; Kil, H. S.; Oh, S. J.; Lee, S. J.; Kim, J. S.; Ryu, J. S.; Moon, D. H.; et al. A New Class of SN2 Reactions Catalyzed by Protic Solvents: Facile Fluorination for Isotopic Labeling of Diagnostic Molecules. J. Am. Chem. Soc. 2006, 128, 16394−16397. (14) Hu, X.; Li, H.; Liang, W.; Han, S. Systematic Study of the Tautomerism of Uracil Induced by Proton Transfer. Exploration of Water Stabilization and Mutagenicity. J. Phys. Chem. B 2005, 109, 5935−5944. (15) Sobolewski, A. L.; Adamowicz, L. Theoretical Investigations of Proton Transfer Reactions in a Hydrogen Bonded Complex of Cytosine with Water. J. Chem. Phys. 1995, 102, 5708−5718. (16) Gorb, L.; Leszczynski, J. Intramolecular Proton Transfer in Monohydrated Tautomers of Cytosine: An Ab Initio Post-Hartree− Fock Study. Int. J. Quantum Chem. 1998, 70, 855−862. (17) Michalkova, A.; Kosenkov, D.; Gorb, L.; Leszczynski, J. Thermodynamics and Kinetics of Intramolecular Water Assisted Proton Transfer in Na+-1-methylcytosine Water Complexes. J. Phys. Chem. B 2008, 112, 8624−8633. (18) Furmanchuk, A.; Isayev, O.; Gorb, L.; Shishkin, O. V.; Hovorun, D. M.; Leszczynski, J. Novel View on the Mechanism of Water-assisted Proton Transfer in the DNA Bases: Bulk Water Hydration. Phys. Chem. Chem. Phys. 2011, 13, 4311−4317. (19) Kim, N. J. DFT Study of Water-assisted Intramolecular Proton Transfer in the Tautomers of Thymine Radical Cation. Bull. Korean Chem. Soc. 2006, 27, 1009−1014. (20) Gu, J.; Leszczynski, J. A DFT Study of the Water-assisted Intramolecular Proton Transfer in the Tautomers of Adenine. J. Phys. Chem. A 1999, 103, 2744−2750. (21) Gorb, L.; Leszczynski, J. Intramolecular Proton Transfer in Mono- and Dihydrated Tautomers of Guanine: An Ab Initio Post Hartree−Fock Study. J. Am. Chem. Soc. 1998, 120, 5024−5032. (22) Gu, J.; Leszczynski, J. Influence of the Oxygen at the C8 Position on the Intramolecular Proton Transfer in C8-Oxidative Guanine. J. Phys. Chem. A 1999, 103, 577−584. (23) Llano, J.; Eriksson, L. A. Oxidation Pathways of Adenine and Guanine in Aqueous Solution from First Principles Electrochemistry. Phys. Chem. Chem. Phys. 2004, 6, 4707−4713. (24) Munk, B. H.; Burrows, C. J.; Schlegel, H. B. Exploration of Mechanisms for the Transformation of 8-Hydroxy Guanine Radical to FAPyG by Density Functional Theory. Chem. Res. Toxicol. 2007, 20, 432−444. (25) Chaban, G. M.; Wang, D.; Huo, W. M. Ab Initio Study of Guanine Damage by Hydroxyl Radical. J. Phys. Chem. A 2015, 119, 377−382. (26) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The Missing Term in Effective Pair Potentials. J. Phys. Chem. 1987, 91, 6269−6271. (27) Valiev, M.; Garrett, B. C.; Tsai, M.-K.; Kowalski, K.; Kathmann, S. M.; Schenter, G. K.; Dupuis, M. Hybrid Approach for Free Energy Calculations with High-level Methods: Application to the SN2 Reaction of CHCl3 and OH− in Water. J. Chem. Phys. 2007, 127, No. 051102.

The calculated free-energy barrier is 19.2 kcal/mol at the DFT/M08-SO/cc-pVTZ+ level of theory. The strongest electrostatic interaction stabilized the proton-transfer process in the transition state; therefore, it lowers the reaction barrier height in aqueous solution. The calculation also shows that the solvent energy contributes −28.5 kcal/mol to the transition state and the polarization effect contributes 19.9 kcal/mol. In total, the water solution contributes −8.6 kcal/mol to the activation barrier, which means that the water solution lowers the transition state in aqueous solution and has a catalytic effect for this proton-transfer reaction mechanism.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.7b09965. Cartesian coordinates (Å) of stationary points (Table S1) (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Dunyou Wang: 0000-0002-5130-3488 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS D.Y. Wang thanks Dr. Chaban for sending us the geometries of the stationary points in gas phase and Dr. Zongliang Li for the advice on spin density. This work was supported by the National Natural Science Foundation of China (Grant Nos. 11774206 and 11374194) and the Taishan Scholarship fund. The computation work was carried out at the Shenzhen Supercomputer Center of China.



REFERENCES

(1) Wang, T.; Yin, H.; Wang, D.; Valiev, M. Hybrid Quantum Mechanical and Molecular Mechanics Study of the SN2 Reaction of CCl4 + OH− in Aqueous Solution: The Potential of Mean Force, Reaction Energetics, and Rate Constants. J. Phys. Chem. A 2012, 116, 2371−2376. (2) Zhang, J.; Xu, Y.; Chen, J.; Wang, D. A Multilayeredrepresentation, Quantum Mechanical/Molecular Mechanics Study of the CH3Cl + F− Reaction in Aqueous Solution: The Reaction Mechanism, Solvent Effects and Potential of Mean Force. Phys. Chem. Chem. Phys. 2014, 16, 7611−7617. (3) Xu, Y.; Zhang, J.; Wang, D. Investigation of the CH3Cl + CN− Reaction in Water: Multilevel Quantum Mechanics/Molecular Mechanics Study. J. Chem. Phys. 2015, 142, No. 244505. (4) Liu, P.; Wang, D.; Xu, Y. A New, Double-inversion Mechanism of the F− + CH3Cl SN2 Reaction in Aqueous Solution. Phys. Chem. Chem. Phys. 2016, 18, 31895−31903. (5) Safi, B.; Choho, K.; Geerlings, P. Quantum Chemical Study of the Thermodynamic and Kinetic Aspects of the SN2 Reaction in Gas Phase and Solution Using a DFT Interpretation. J. Phys. Chem. A 2001, 105, 591−601. (6) Bathgate, R. H.; Moelwyn-Hughes, E. A. The Kinetics of Certain Ionic Exchange Reactions of the Four Methyl Halides in Aqueous Solution. J. Chem. Soc. 1959, 2642−2648. (7) Fu, A.-p.; Li, H.-l.; Du, D.-m.; Zhou, Z.-y. Theoretical Study on the Reaction Mechanism of Proton Transfer in Formamide. Chem. Phys. Lett. 2003, 382, 332−337. 3131

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132

Article

The Journal of Physical Chemistry B (28) Yin, H.; Wang, D.; Valiev, M. Hybrid Quantum Mechanical/ Molecular Mechanics Study of the SN2 Reaction of CH3Cl + OH− in Water. J. Phys. Chem. A 2011, 115, 12047−12052. (29) Zheng, J.; Zhao, Y.; Truhlar, D. G. The DBH24/08 Database and Its Use to Assess Electronic Structure Model Chemistries for Chemical Reaction Barrier Heights. J. Chem. Theory Comput. 2009, 5, 808−821. (30) Fox, T.; Kollman, P. A. Application of the RESP Methodology in the Parametrization of Organic Solvents. J. Phys. Chem. B 1998, 102, 8070−8079. (31) Valiev, M.; Bylaska, E. J.; Govind, N.; Kowalski, K.; Straatsma, T. P.; Van Dam, H. J. J.; Wang, D.; Nieplocha, J.; Apra, E.; Windus, T. L.; et al. NWChem: A Comprehensive and Scalable Open-source Solution for Large Scale Molecular Simulations. Comput. Phys. Commun. 2010, 181, 1477−1489. (32) Zhang, Y.; Liu, H.; Yang, W. Free Energy Calculation on Enzyme Reactions with an Efficient Iterative Procedure to Determine Minimum Energy Paths on a Combined Ab Initio QM/MM Potential Energy Surface. J. Chem. Phys. 2000, 112, 3483−3492. (33) Henkelman, G.; Uberuaga, B. P.; Jónsson, H. A Climbing Image Nudged Elastic Band Method for Finding Saddle Points and Minimum Energy Paths. J. Chem. Phys. 2000, 113, 9901−9904. (34) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A. et al. Gaussian 09, revision D.01; Gaussian, Inc.: Wallingford, CT, 2013. (35) Hu, H.; Lu, Z.; Parks, J. M.; Burger, S. K.; Yang, W. Quantum Mechanics/Molecular Mechanics Minimum Free-energy Path for Accurate Reaction Energetics in Solution and Enzymes: Sequential Sampling and Optimization on the Potential of Mean Force Surface. J. Chem. Phys. 2008, 128, No. 034105. (36) Liu, P.; Wang, Q.; Niu, M.; Wang, D. Multi-level Quantum Mechanics and Molecular Mechanics Study of Ring Opening Process of Guanine Damage by Hydroxyl Radical in Aqueous Solution. Sci. Rep. 2017, 7, No. 7798.

3132

DOI: 10.1021/acs.jpcb.7b09965 J. Phys. Chem. B 2018, 122, 3124−3132