PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
J Am Chem Soc. Author manuscript; available in PMC 2010 May 27.
Published in final edited form as:
PMCID: PMC2710000
NIHMSID: NIHMS108683

Theoretical Analysis of the Unusual Temperature Dependence of the Kinetic Isotope Effect in Quinol Oxidation

Abstract

This paper presents theoretical calculations on model biomimetic systems for quinol oxidation. In these model systems, an excited-state [Ru(bpy)2(pbim)]+ complex (bpy = 2,2’-dipyridyl, pbim = 2-(2-pyridyl)benzimidazolate) oxidizes a ubiquinol or plastoquinol analogue in acetonitrile. The charge transfer reaction occurs via a proton-coupled electron transfer (PCET) mechanism, in which an electron is transferred from the quinol to the Ru and a proton is transferred from the quinol to the pbim ligand. The experimentally measured average kinetic isotope effects (KIEs) at 296 K are 1.87 and 3.45 for the ubiquinol and plastoquinol analogues, respectively, and the KIE decreases with temperature for plastoquinol but increases with temperature for ubiquinol. The present calculations provide a possible explanation for the differences in magnitudes and temperature dependences of the KIEs for the two systems and, in particular, an explanation for the unusual inverse temperature dependence of the KIE for the ubiquinol analogue. These calculations are based on a general theoretical formulation for PCET reactions that includes quantum mechanical effects of the electrons and transferring proton, as well as the solvent reorganization and proton donor-acceptor motion. The physical properties of the system that enable the inverse temperature dependence of the KIE are a stiff hydrogen bond, which corresponds to a high-frequency proton donor-acceptor motion, and small inner-sphere and solvent reorganization energies. The inverse temperature dependence of the KIE may be observed if the (0/0) pair of reactant/product vibronic states is in the inverted Marcus region, while the (0/1) pair of reactant/product vibronic states is in the normal Marcus region and is the dominant contributor to the overall rate. In this case, the free energy barrier for the dominant transition is lower for deuterium than for hydrogen because of the smaller splittings between the vibronic energy levels for deuterium, and the KIE increases with increasing temperature. The temperature dependence of the KIE is found to be very sensitive to the interplay among the driving force, the reorganization energy, and the vibronic coupling in this regime.

I. Introduction

Redox reactions of quinol species play an important role in a wide range of biological energy transduction processes.1 For example, the oxidation of ubiquinol by the cytochrome bc1 complex is a key step in the catalytic cycle that transports electrons while pumping protons across a membrane in respiration and photosynthesis.2 The cyt bc1 complex may also be used for superoxide production in mitochondria.3 Understanding the mechanism of ubiquinol oxidation by the cyt bc1 complex and related systems is important for both biomedical and energy-based technological applications.

Recently Kramer and coworkers studied the kinetics of quinol oxidation by the cyt bc1 complex and biomimetic model systems in acetonitrile.3 For both the biological and biomimetic systems, they found that the oxidation of ubiquinol exhibits an unusual temperature dependence of the kinetic isotope effect (KIE), where the KIE increases with increasing temperature (i.e., the apparent activation free energy is greater for hydrogen than for deuterium). This temperature dependence is considered to be unusual because most semiclassical models based on transition state theory4, 5 predict that both zero point energy and hydrogen tunneling effects result in a smaller activation free energy for hydrogen than for deuterium. For both the biological and biomimetic systems, they found that plastoquinol oxidation does not exhibit this unusual temperature dependence of the KIE. Similar behavior was observed by Nagaoka and coworkers for the oxidation of ubiquinol by a tocopherol radical in ethanol.6, 7 Thus, the unusual inverse temperature dependence of the KIE has been observed for quinol oxidation in three different reaction systems.

In this paper, we perform theoretical calculations on the model biomimetic complexes studied experimentally by Kramer and coworkers.3 The oxidant is an excited-state [Ru(bpy)2(pbim)]+ complex (bpy = 2,2’-dipyridyl, pbim = 2-(2-pyridyl)benzimidazolate) in acetonitrile. The ubiquinol and plastoquinol analogues are 2,3-dimethoxy-5-methyl-1,4-benzoquinol (UQH2) and 2,3,5-trimethyl-1,4-benzoquinol (PQH2), respectively. In the experiments, the Ru complex is photoexcited to a metal-to-ligand charge transfer (MLCT) state, and the emission quenching kinetics of the excited-state complex are analyzed. The data suggest that the quinol is hydrogen-bonded to the pbim ligand of the Ru complex, and the charge transfer reaction occurs via a proton-coupled electron transfer (PCET) mechanism810 in which an electron is transferred from the quinol to the Ru(III) and a proton is transferred from the quinol to the pbim ligand. This mechanism is illustrated in Figure 1, and the UQH2 and PQH2 molecules are depicted in Figure 2. The electron and proton are expected to transfer simultaneously to avoid the high-energy intermediates involved in the sequential reactions. Previously, the single proton and single electron transfer reaction were estimated to be endothermic by ~40 kJ/mol.3 Further evidence of a PCET mechanism is that the average KIEs at 296 K are 1.87 and 3.45 for the UQH2 and PQH2 systems, respectively.3 As mentioned above, the KIE increases with temperature for the UQH2 system and decreases with temperature for the PQH2 system. Yamamoto and Kato studied a simplified model of the UQH2 oxidation reaction theoretically with quantum chemical methods.11 While their studies suggest the significance of a stiffer proton potential at the transition state, they were unable to reproduce or explain the experimentally observed unusual inverse temperature dependence of the KIE for the UQH2 system.

Figure 1
PCET reaction corresponding to oxidation of UQH2 following photoexcitation of the [Ru(bpy)2(pbim)]+ complex to a MLCT state. The proton transfer and electron transfer reactions are indicated with arrows. The overall charge of the hydrogen-bonded complex ...
Figure 2
(a) UQH2 and (b) PQH2.

The calculations in the present paper are based on a general theoretical formulation for PCET reactions that includes quantum mechanical effects of the electrons and transferring proton, as well as the solvent reorganization and proton donor-acceptor motion.1215 Alternative theoretical treatments of PCET and vibrationally nonadiabatic proton transfer reactions have been used to investigate the temperature dependence of KIEs5, 1623 but, to our knowledge, none of these studies has revealed an inverse temperature dependence, where the KIE increases with increasing temperature. Our calculations provide a possible explanation for the differences in magnitudes and temperature dependences of the KIEs for the UQH2 and PQH2 systems and, in particular, an explanation for the inverse temperature dependence of the KIE for the UQH2 system. The physical characteristics of the system that enable the inverse temperature dependence are a stiff hydrogen bond, which corresponds to a high-frequency proton donor-acceptor mode, and small inner-sphere and solvent reorganization energies. The unusual temperature dependence can arise when the (0/0) pair of reactant/product vibronic states is in the inverted Marcus region and excited vibronic states contribute significantly to the overall reaction rates. The temperature dependence of the KIE is found to be very sensitive to the interplay among the driving force, the reorganization energy, and the vibronic coupling.

An outline of the paper is as follows. In Section II, we summarize the theoretical formulation for PCET and the computational methods used to calculate the input quantities for the rate constant expressions. In Section III, we present the results and analysis. The conclusions are presented in Section IV.

II. Theory and Methods

A. Rate constant expressions

A series of nonadiabatic rate constant expressions for PCET reactions have been derived previously.1215 In this framework, the PCET reaction is described in terms of nonadiabatic transitions between the reactant and product electron-proton vibronic states. The derivations and discussions of the underlying approximations and regimes of validity are given elsewhere.1215 Here we discuss only the rate constant expressions that are directly relevant to this system. For fixed proton donor-acceptor distance R, the PCET rate constant is13

k=μPμν|Vμν|2πλμνkBTexp[(ΔGμν0+λμν)24λμνkBT]
(1)

where the summations are over reactant and product vibronic states, Pμ is the Boltzmann probability for the reactant state μ Vµν is the vibronic coupling between the reactant and product vibronic states μ and ν, λµν is the reorganization energy for states μ and ν, and ΔGμν0 is the free energy of reaction for states μ and ν. The quantities λµν and ΔGμν0 are defined analogously to the traditional definitions in Marcus theory for electron transfer24 and can be calculated from the vibronic free energy surfaces corresponding to states μ and ν using the explicit definitions given in Ref 13. As discussed elsewhere,15 the vibronic coupling is often approximated as Vµν≈ νelSµν, where νel is the electronic coupling and Sµν is the overlap integral between the reactant and product proton vibrational wavefunctions. This approximation is valid in the limit of electronically nonadiabatic proton transfer.25, 26 The rate constant expression in Eq. (1) is expected to be applicable to systems with very stiff hydrogen bonds for which the motion along the R coordinate is restricted.

Analytical rate constant expressions including the dynamical effects of the R motion have also been derived.14, 15 In these derivations, the vibronic coupling between pairs of states is assumed to depend exponentially on R: Vμν=Vμν(0)exp[αμν(RR¯)], where R is the equilibrium value of R assumed to be the same for all states, Vμν(0) is the vibronic coupling between states μ and ν at a distance R, and αµν is a parameter that reflects the distance dependence of the vibronic coupling. In the limit of electronically nonadiabatic proton transfer, Vμν(0)=VelSμν(0), where Sμν(0) is the overlap at the distance R and the parameter αµν describes the approximately exponential decay of the overlap with R in the region near R = R. In this limit, the parameters αµν are determined for each pair of vibronic states by calculating the numerical derivative of the natural logarithm of the corresponding overlap integral with respect to R at R = R.

Using the short-time, high-temperature approximation for the solvent modes and representing the R-mode time correlation function by that of a quantum mechanical harmonic oscillator, the rate constant can be expressed as:14

k=μPμν|Vμν(0)|22Ωexp[2λμν(α)ζΩ]dτexp[12χτ2+p(cosτ1)+i(qsinτ+θτ)]
(2)

with the dimensionless parameters defined as

ζ=coth(12βΩ); χ=2λβ2Ω2; θ=ΔGμν0+λΩ p=ζλμν(α)Ω; q=λμν(α)Ω
(3)

Here β=1/kBT, λμν(α) is defined as λμν(α)=2αμν2/2M, and M and Ω are the R-mode effective mass and frequency, respectively. Note that the reorganization energy λ is assumed to be the same for all pairs of states in this expression. The free energy of reaction is often expressed as ΔGμν0=ΔG0+Δεμν, where ΔG0ΔG000 and Δεµν is the difference between the product and reactant vibronic energy levels relative to their respective ground vibronic states.

Simplified analytical rate constant expressions have been derived in the low-frequency and high-frequency limits for the R-mode.14 As will be discussed below, density functional theory (DFT) calculations indicate that the quinol systems studied in this paper have a relatively high-frequency R-mode. Moreover, in the low-frequency R-mode limit, where ħΩ<<kBT the KIE has the approximate form27 kHkDexp[2kBTMΩ2(αD2αH2)]|SH|2|SD|2 when only the lowest energy reactant and product vibronic states (i.e., μ = ν = 0) are included, leading to a decrease in the KIE with increasing temperature since αD > αH.15 Even when excited vibronic states are included, the exponential temperature-dependent term tends to dominate the temperature dependence of the KIE, and the KIE decreases with increasing temperature. Thus, the unusual inverse temperature dependence of the KIE observed for the UQH2 system is not predicted by the low-frequency R-mode rate constant expression. This analysis of the form of the KIE in the low-frequency R-mode limit, in conjunction with the DFT prediction of a high-frequency R-mode, imply that these systems are in the high-frequency R-mode regime. However, the derivation of the simplified expression in the high-frequency R-mode limit is based on the assumption that λ>|ΔGμν0| for all relevant pairs of reactant/product vibronic states,14 and this condition does not hold for these systems. Thus, we use the full expression given in Eq. (2), henceforth denoted the dynamical rate constant expression, and show that it leads to similar results as the fixed-R rate constant expression given in Eq. (1) for these systems.

The rate constant expressions in Eqs. (1) and (2) correspond to unimolecular PCET within the hydrogen-bonded complex shown in Figure 1, whereas the rate constants measured experimentally for these systems were second-order bimolecular rate constants. The relation between the bimolecular and unimolecular rate constants is kbi=Keq kuni, where Keq is the equilibrium association constant to form the hydrogen-bonded complex. In our analysis, we assume that Keq is independent of isotope. This approximation is justified by previous DFT calculations of equilibrium isotope effects. As given in the Supporting Information of Ref 3, the equilibrium isotope effect for intramolecular hydrogen bonding in UQH2 was calculated to be ~1.01. Moreover, the equilibrium isotope effects for hydrogen bond formation in ubiquinol-imidazolate complexes were also calculated to be negligible.28

B. Calculation of input quantities

For these PCET systems, the four diabatic states corresponding to the relevant charge transfer states are defined as:12

QOHNRu(III)
(1a)

QOHNRu(III)
(1b)

QOH·+NRu(II)
(2a)

QO·HNRu(II)
(2b)

where 1 and 2 denote the electron transfer state, and a and b denote the proton transfer state. We define the reactant electronic state to be a mixture of the (1a) and (1b) diabatic states and the product electronic state to be a mixture of the (2a) and (2b) diabatic states. In our implementation, the proton vibrational states are calculated for the proton potential energy curves associated with the reactant and product electronic states, and the reactant and product vibronic states are products of the corresponding electronic and proton vibrational states.

Various input quantities are required to determine the rates and KIEs. The atomic coordinates and the partial charges for the four diabatic states are required for the calculation of the solvent (outer-sphere) reorganization energies. The equilibrium proton donor-acceptor distance and the associated effective frequency are also input quantities for the rate constant expressions. In addition, the proton potential energy curves are required for the calculation of the proton vibrational wavefunction overlaps in the vibronic couplings. We performed a series of electronic structure calculations on model systems to obtain these input quantities. All electronic structure calculations discussed in this paper were performed with Gaussian03.29

The objective of this study is not to obtain quantitatively accurate values for these input quantities, but rather to obtain qualitatively reasonable values that enable us to model the experimental data. Thus, we perform relatively low-level calculations to obtain estimates of the solvent reorganization energies, driving forces, and proton donor-acceptor frequencies. These estimates provide an indication of the range of parameters that should be considered. Subsequently, since the magnitude and temperature dependence of the KIEs are very sensitive to these values, we treat these quantities as parameters to fit the experimental data. The parameter values resulting from the fitting procedure should still be viewed as qualitative estimates because they depend strongly on the proton potentials, which are also only qualitatively accurate for these systems. Despite these limitations in terms of quantitative accuracy, this modeling procedure provides a plausible explanation for the experimentally observed inverse temperature dependence of the KIE for UQH2 and for the differences between the KIEs for UQH2 and PQH2.

1. Atomic coordinates and charges

We calculated the atomic coordinates using three different models. For the full system gas phase model, we optimized the geometries for the Ru(bpy)2(pbim)+ complex hydrogen bonded to either UQH2 or PQH2 in the gas phase with density functional theory (DFT) at the B3LYP/LANL2DZ level.3033 For the truncated models, we optimized the geometries for the pbim ligand hydrogen bonded to either UQH2 or PQH2 with DFT at the B3LYP/6-31G** level in the gas phase and also using the polarized continuum model (PCM)34, 35 to represent the acetonitrile solvent. These two models are denoted the truncated gas phase and truncated solvated models, respectively. A comparison of the O–N distances at the hydrogen-bonding interface for these three models is given in Table 1. Note that this distance is nearly identical for the full system gas phase and truncated solvated models. We used the O–N distance from the truncated solvated model as the equilibrium proton donor-acceptor distance in our PCET rate constant calculations. The structure used to calculate the solvent reorganization energies and the proton potential energy curves was obtained by superimposing the truncated solvated model, which consists of the quinol and pbim ligand in acetonitrile, onto the full system gas phase model. Note that the structure of the Ru complex was optimized for the ground electronic state, although the reaction occurs in the MLCT excited state. The resulting differences in the coordinates are not expected to significantly impact the calculation of the solvent reorganization energies, which depend mainly on the change in solute charge distribution. As discussed below, the partial atomic charges are chosen to represent the MLCT excited state. Moreover, these calculated solvent reorganization energies are viewed as qualitative estimates, and subsequently this quantity is treated as a parameter that is fit to the experimental data.

Table 1
The O–N equilibrium proton donor–acceptor distance in Å determined from DFT calculations for UQH2 and PQH2 with the three different models defined in the text.

We also calculated the normal mode frequencies for the full system gas phase model and identified the mode dominated by the proton donor-acceptor (O–N) motion. The frequency Ω and reduced mass M of this mode were 1003 cm−1 and 7.22 amu for UQH2 and 1000 cm−1 and 7.02 for PQH2. Note that the reduced masses are similar to the reduced mass of 7.47 amu for oxygen and nitrogen. The small mass and high frequency indicate that this mode is localized to the proton transfer interface. Moreover, the magnitudes of these frequencies indicate that these systems are in the high-frequency R-mode regime, which is defined as ħΩ>>kBT. In addition, these frequencies are very similar for the two systems and hence are probably not responsible for differences in the magnitudes and temperature dependences of the KIEs.

The atomic charges for the four diabatic states given above are required for the calculation of the solvent reorganization energies. We determined these atomic charges by performing separate calculations on the isolated ligands and quinols with overall charges corresponding to the diabatic states. For simplification, we placed an electron on one of the bpy ligands to mimic the MLCT state, although the electron may be delocalized over both bpy ligands. This approximation does not significantly impact the solvent reorganization energies because the atomic charges on the bpy ligands are the same for all four diabatic states. We obtained the atomic charges for the protonated and deprotonated quinols, the protonated pbimH ligand, the deprotonated pbim ligand, and the bpy and bpy ligands by optimizing the isolated ligands with DFT at the B3LYP/6-31G** level and applying the CHELPG method36 to the optimized ligands. The ruthenium atom was assigned a charge of +3 or +2 corresponding to the appropriate oxidation state. Note that the application of the CHELPG method to the full complex in a manner that is consistent for all four diabatic states, including those corresponding to electronically excited states, is not practical. To maintain the proper charge relations among the four diabatic states,12 the atomic charges for the 2a diabatic state were obtained as a function of the atomic charges in the other three diabatic states (i.e., q2a = q2bq1b + q1a).

2. Reorganization energies and driving forces

The solvent (outer-sphere) reorganization energies were calculated with the frequency-resolved cavity model (FRCM)37, 38 within the framework of the multistate continuum theory for PCET.12 This approach allows for distinct effective solute cavities pertaining to the optical and inertial solvent response. The cavities are formed from spheres centered on all of the atoms. The two effective radii for the solute atoms are defined as rrvdW and rin =r+δ, where rvdW is the van der Waals radius, κ is a universal scaling factor, and δ is a parameter related to the size of the solvent molecules. The factor for the van der Waals radius and the thickness of the inner layer were κ = 0.9 and δ = 1.8, respectively,38 and the static and optical dielectric constants were εo = 37.5 and ε = 1.8, respectively, for acetonitrile at 296 K. The temperature dependence of the dielectric constants is neglected because this effect is relatively insignificant for the temperature range studied.

Using methodology described previously,12 the solvent reorganization energies for the (0/0) pair of reactant/product vibronic states were calculated to be 6.35 kcal/mol and 6.29 kcal/mol for the UQH2 and PQH2 systems, respectively. These values are relatively low mainly because the reaction from state 1a to 2b corresponds to a net electron transfer from the nitrogen atom of the pbim ligand to the Ru, corresponding to an electron transfer over 4.2 Å along a pathway buried in the solute cavity screened from solvent. In addition, acetonitrile has a lower dielectric constant than water and a larger δ parameter, which dictates the size of the outer cavity, due to the larger molecular size of acetonitrile compared to water. Moreover, these estimates should be viewed as upper limits for the solvent reorganization energies because the choice of atomic charges (i.e., the assignment of +3 or +2 to the ruthenium atom) neglects charge transfer between the ruthenium and the ligands. As shown previously,39 this simplification to the charge distribution does not qualitatively alter the calculated solvent reorganization energies but could lead to quantitative overestimation of the solvent reorganization energies.

The inner-sphere reorganization energy corresponding to the solute modes can be estimated using standard approaches.24, 39 Experimental results, however, imply that the inner-sphere reorganization energy does not contribute significantly to these types of reactions. In particular, crystal structures of Ru(bpy)33+ and Ru(bpy)32+ indicate that the Ru–N distances are the same to within experimental error.40 Moreover, the inner-sphere reorganization energy within the quinol is also expected to be relatively small.

The driving forces for PCET from UQH2 and PQH2 to the Ru complex were estimated to be −6.0 kcal/mol and −4.5 kcal/mol, respectively, using experimental pKa values and redox potentials. The physical basis for these estimates is given in Supporting Information, but they should be viewed as only qualitative approximations. Due to the uncertainties in the estimates of reorganization energies and driving forces, we treat these quantities as parameters that are fit to the experimental data. Nevertheless, these estimates provide an indication of the ranges of parameters that are physically meaningful for these systems. Moreover, the similarities in the magnitudes of the reorganization energies and driving forces suggest that the system is in the vicinity of the inverted region, where –ΔGo > λ.

3. Proton potentials

We generated the gas phase proton potential energy curves by moving the hydrogen nucleus on a grid along the O–N axis and calculating the energy for the grid points with DFT at the B3LYP/6-31G** level. As discussed above, the reactant and product states are defined for the complex after photoexcitation to the MLCT state, so the ruthenium is in the +3 oxidation state for the reactant and the +2 oxidation state for the product. In order to calculate the reactant and product proton potential energy curves consistently, the energies were calculated for a reduced model in which the Ru and bpy ligands were removed from the superimposed full gas phase structure, and a point charge was placed at the Ru site to mimic the charge of the Ru and the negatively charged bpy ligand following photoexcitation. For the reactant structure, the quinol-pbim subsystem was assigned a charge of −1, and a point charge of +2 was placed at the Ru site (i.e., the sum of +3 for the Ru and −1 for the bpy ligand). For the product structure, the quinol-pbim subsystem was assigned a charge of 0, and a point charge of +1 was placed at the Ru site (i.e., the sum of +2 for the Ru and −1 for the bpy ligand). All nuclei were fixed for the calculation of the proton potential energy curves except the transferring hydrogen nucleus.

We calculated the solvated proton potential curves using the multistate continuum theory12 in conjunction with the FRCM method.37, 38 First we fit the gas phase proton potential energy curves with the molecular mechanical functional form used in Ref. 41. The parameters are given in Supporting Information. Then we calculated the lowest energy reactant and product state vibronic free energy surfaces as functions of two collective solvent coordinates corresponding to proton and electron transfer. Subsequently, we generated the solvated reactant and product proton potentials at the collective solvent coordinates corresponding to the intersection of these vibronic surfaces. The resulting solvated proton potentials were represented numerically by a spline fit, and the proton vibrational wavefunctions corresponding to these solvated proton potentials were calculated numerically with Fourier grid methods. The gas phase and solvated proton potentials are depicted in Figure 3 for both the ubiquinol and plastoquinol systems. Note that the proton potentials for UQH2 and PQH2 are virtually indistinguishable for the product state but are noticeably different for the reactant state.

Figure 3
Proton potential energy curves for UQH2 (left) and PQH2 (right). Reactant curves are shown on the top in blue, and product curves are shown on the bottom in red. Gas phase curves are dashed, and solvated curves are solid.

III. Results and Discussion

As mentioned above, the experimentally observed KIE is 1.87 for UQH2 and 3.45 for PQH2 at 296 K.3 Table 1 indicates that the proton donor-acceptor distance is 2.65 Å for UQH2 and 2.67 Å for PQH2. Our previous analyses of model systems suggest that the larger KIE for PQH2 is due mainly to the slightly larger proton donor-acceptor distance for PQH2.15 Typically the KIE increases as the proton donor-acceptor distance increases for fixed other parameters because the KIE is approximately proportional to the square of the ratio of the hydrogen and deuterium vibrational overlaps for a given pair of states (i.e., KIESH2/SD2). This ratio increases with increasing proton donor-acceptor distance because the deuterium overlap decreases faster than the hydrogen overlap as this distance increases.

Furthermore, as discussed above, the experimentally measured KIE increases with temperature for UQH2 and decreases with temperature for PQH2.3 To reproduce the experimentally observed KIEs for both quinol systems, we implemented a parameter fitting procedure using both the fixed-R rate constant expression and the expression including the dynamical effects of the R motion. In this procedure, we fit the KIEs to the experimental data by minimizing the sum of the squares of the differences between the experimental and calculated KIEs. We varied λ and ΔG0 for the fixed-R rate constant expression in Eq. (1), and we varied λΔG0, and Ω for the dynamical rate constant expression in Eq. (2). For the latter calculations, we used the reduced mass of oxygen and nitrogen, M = 7.47 amu, for the effective mass associated with the R-mode. For simplicity, we used the same value of ΔG0ΔG000 for both H and D for each system. In principle, this value should differ by

ΔΔZPE=[ZPEDproductZPEDreactant][ZPEHproductZPEHreactant]
(4)

, where ZPEHreactant is defined as the zero point energy for hydrogen moving in the reactant proton potential and the other quantities are defined analogously. We calculated ΔΔZPE to be ~0.1 kcal/mol for both systems, indicating that this contribution is negligible for our studies.

The parameter values are given in Table 2. For both systems, the reorganization energy decreases relative to the value calculated with the FRCM method. This difference is expected because the partial atomic charges assigned to the four diabatic states used in the FRCM calculations are too localized, leading to an overestimation of the reorganization energy. The magnitude of the driving force also decreases for both systems. This difference is also expected because the estimates of the driving force were based on pKa values and redox potentials and hence assumed infinite separation of the reactant and product species, leading to an overestimation of the driving forces. However, the trend in the driving forces is still retained, namely that the magnitude of the driving force is greater for the UQH2 system than for the PQH2 system. The frequencies are significantly larger than the frequencies calculated from the DFT gas phase calculations. This difference may be due in part to the neglect of solvent effects in the DFT frequency calculations.

Table 2
Parameters obtained from fitting procedure described in the text. The fixed-R rate constant expression is given in Eq. (1), and the dynamical rate constant expression is given in Eq. (2). The mass is M = 7.47 amu for the calculations using Eq. (2). The ...

We emphasize that these fitted values of the parameters are not quantitatively meaningful because of the simplicity of the model. Due to the complex interplay between the vibronic coupling and these parameters, the quantitative values of these parameters depend on the proton potential energy curves, which were generated with an approximate method. Thus, the unusually low reorganization energies and the unusually high proton donor-acceptor frequencies are not quantitatively meaningful and simply imply a qualitative trend. The analysis of the temperature dependence of the KIE given below does not rely on quantitatively accurate parameters but rather depends on only these qualitative trends.

As shown in Figure 4, we are able to reproduce the experimentally observed magnitude and temperature dependence of the KIE for both quinols using this fitting procedure. This figure was generated with the fixed-R rate constant expression in Eq. (1), but very similar results were obtained using the dynamical rate constant expression in Eq. (2). Figure 4 illustrates that the KIE increases with temperature for UQH2 and decreases with temperature for PQH2. Both systems have a small solvent reorganization energy and a high-frequency R-mode, which corresponds to a stiff hydrogen bond. For both reactions, the (0/0) pair of reactant/product vibronic states is in the inverted Marcus region, where -ΔG0 > λ. Under these conditions, we have found that the KIE could either increase or decrease with temperature, depending on a subtle interplay among the various parameters. For these two quinol systems, the qualitatively different temperature dependences of the KIE are due mainly to the differences in the driving forces and proton potentials.

Figure 4
Plot of ln[KIE} versus 1/T calculated with the fixed-R rate constant expression (open circles) using the parameters in Table 2 and the experimental values (closed circles) for (a) UQH2 and (b) PQH2.

To elucidate the differences between these two systems, we analyzed the contributions of the pairs of reactant/product vibronic states using the fixed-R rate constant expression given in Eq. (1). The calculated free energy curves for the lowest three reactant and product vibronic states, along with the associated proton vibrational wavefunctions, are depicted in Figure 5. The relative contribution of each pair of reactant/product vibronic states is determined by a balance among the Boltzmann probability Pμ for the reactant state, the square of the vibronic coupling Vμν2Sμν2, and the exponential factor exp [ΔGμν/kBT], where the free energy barrier can be calculated as ΔGμν=(ΔGμν0+λμν)2/4λμν. The balance among these factors, and therefore the relative contributions from each pair of reactant/product vibronic states, depends on temperature and differs for H and D, leading to complex behavior of the temperature dependence of the KIE. We examined the contributions from the pairs of reactant/product vibronic states for UQH2 and PQH2 with H and D at two different temperatures. The analysis at 296 K is provided in Table 3, and the analysis at 320 K is provided in Supporting Information.

Figure 5
Free energy curves as functions of a collective solvent coordinate (center) with the proton potential energy curves and hydrogen vibrational wavefunctions on the left (reactant) and on the right (product) for (a) UQH2 and (b) PQH2. The three lowest energy ...
Table 3
Analysis of the contributions of pairs of reactant/product vibronic states µ/ν for UQH2 and PQH2 systems at 296 K with parameters given in Table 2. The % Cont. is the percentage contribution from the specified pair of reactant/product ...

For both systems, the (0/1) pair of reactant/product vibronic states contributes substantially to the overall rate. The free energy barrier exponential term favors the (0/0) pair, and the vibronic coupling factor favors the (0/1) pair. As shown in Figure 5, the (0/0) pair is in the inverted region, and the (0/1) pair is in the normal region. As a result, the free energy barriers for these two transitions are similar, with the (0/1) free energy barrier only slightly higher than the (0/0) free energy barrier. The vibronic coupling factor favors the (0/1) pair because the overlap between the reactant and product proton vibrational wavefunctions is greater for the excited product vibronic state. Due to the similar free energy barriers for the (0/0) and (0/1) pairs, their relative contributions to the overall rate are significantly influenced by the vibronic couplings, which favor the (0/1) transition.

The analysis in Table 3 indicates a significant qualitative difference between the ubiquinol and plastoquinol systems. For hydrogen transfer, UQH2 has a dominant contribution from the (0/1) pair of reactant/product vibronic states, while PQH2 has more similar contributions from the (0/0) and (0/1) pair of reactant/product vibronic states. For UQH2, the (0/0) pair of reactant/product states contributes 10.6% and the (0/1) pair contributes 88.9% to the overall rate, while for PQH2, the (0/0) pair contributes 39.1% and the (0/1) pair contributes 59.6% to the overall rate. As discussed above, the free energy barrier exponential term favors the (0/0) pair, and the vibronic coupling factor favors the (0/1) pair. For the PQH2 system, however, the penalty from the increased free energy barrier for the (0/1) pair relative to the (0/0) pair is greater because the increase in free energy barrier is 2.69 kcal/mol for PQH2 and 1.10 kcal/mol for UQH2. This difference arises from relatively subtle differences in the free energies of reaction (i.e., the driving forces) and the proton potential energy curves.

Table 3 also indicates that the (1/2) pair of reactant/product vibronic states contributes significantly for deuterium transfer in both systems. This pair contributes more than other pairs involving the first excited reactant vibronic state because the free energy barrier is nearly zero for this pair. For hydrogen transfer, the contributions from this pair are nearly negligible because of the small Boltzmann probability for the reactant state. The (1/2) pair contributes more for deuterium than for hydrogen because the vibronic energy level splittings are smaller for deuterium, leading to a greater Boltzmann probability for the first excited reactant state. For deuterium, the contribution from the (1/2) pair is greater for PQH2 than for UQH2 because the relative vibronic couplings (i.e., vibrational wavefunction overlaps) favor the (1/2) pair more for PQH2, mainly due to subtle differences in the proton potential energy curves. For both systems, increasing the temperature increases the contributions from the (1/2) pair mainly because of the increased Boltzmann probability for the first excited reactant state.

The unusual temperature dependence for UQH2 can be understood by analyzing the terms in Table 3. This analysis is based on the rate constant expressions valid in the high-frequency R-mode regime. Since the (0/1) pair is dominant for H and D at both 296 and 320 K, we focus on only this pair in the analysis. For this pair of states, the KIE is approximately proportional to the ratio of exp exp[ΔG01/kBT]S012 for H and D, assuming that the reorganization energies are independent of isotope. The vibronic couplings, which are proportional to the square of the overlaps given in this table, are assumed to be independent of temperature. In this case, the temperature dependence of the KIE is determined mainly by the factor exp [(ΔG01(H)ΔG01(D))/kBT], which will increase as the temperature increases when ΔG01(H)>ΔG01(D). As illustrated in Figure 6, the (0/0) pair of reactant/product vibronic states is in the inverted Marcus region, while the (0/1) pair of reactant/product vibronic state is in the normal Marcus region. The reaction free energy ΔG000 is virtually identical for H and D because the differences in zero point energies for the reactant and product cancel (i.e., ΔΔZPE is negligible). In this case, the reaction free energy ΔG010, and therefore the free energy barrier ΔG01, is lower for D than for H because the vibronic energy level splittings are smaller for D. Specifically, ΔG01(H)=2.89 and ΔG01(D)=1.01. Since ΔG01(H)>ΔG01(D), and the (0/1) pair is the dominant contributor to the overall rate, the KIE increases as the temperature increases. Thus, the unusual inverse temperature dependence is easily explained within this theoretical framework.

Figure 6
Free energy curves as functions of a collective solvent coordinate for the UQH2 system. The lowest energy reactant vibronic state (blue, labeled μ=0) and lowest two product vibronic states (red, labeled ν=0 and ν=1) are shown for ...

The analysis of the temperature dependence of the KIE for PQH2 is more complicated because of the significant contributions from more pairs of states. As a result, the KIE is not simply the ratio of factors for a single pair of states, but rather is a ratio of sums over pairs of states. Moreover, since the (1/2) pair contributes significantly for D, the temperature dependence of the KIE is influenced by the Boltzmann probability Pμ as well as the factor exp [ΔGμν/kBT]. An analysis of only the (0/1) pair for the PQH2 system would predict an increase in the KIE as temperature increases for the reasons described above for UQH2. For the PQH2 system, however, the significant contributions from the (0/0) pair for H and the (1/2) pair for D lead to a decrease in the KIE with increasing temperature. In other words, the apparent activation free energies for H and D arise from a linear combination of two or more exponential terms. These calculations imply that the greater electron donating character of the OCH3 substituents in UQH2 relative to the CH3 substituents in PQH2 impacts the driving forces and the proton potential energy curves enough to significantly alter the relative contributions from the various pairs of vibronic states, thereby leading to different temperature dependences of the KIE.

We emphasize that the objective of these calculations was to provide an explanation for the unusual inverse temperature dependence of the KIE for ubiquinol oxidation. In our modeling procedure, the reorganization energies and driving forces are treated as empirical parameters, and small modifications of these parameters may compensate for approximations invoked in other parts of the modeling, such as the generation of the proton potential energy curves. Although the specific values of these parameters are not expected to be quantitatively accurate, the qualitative picture emerging from these studies provides an explanation for the unusual inverse temperature dependence of the KIE. Thus, despite the many underlying approximations, this modeling procedure enabled us to uncover the fundamental mechanism of quinol oxidation in these systems and to assist in the interpretation of the experimental data.

An alternative explanation for the experimental observations is that the precursor complex formation could depend on both temperature and isotope in a manner that leads to this unusual temperature dependence of the KIE. As discussed above, however, the equilibrium isotope effects for hydrogen bond formation in related model systems and for intramolecular hydrogen bonding in UQH2 were previously calculated to be negligible. Note that our explanation of these experimental observations relies on a high-frequency proton donor-acceptor motion. This type of inverse temperature dependence of the KIE is not expected to be observed when the proton donor-acceptor frequency is low relative to the thermal energy. In this low-frequency regime, the temperature dependence of the KIE is dominated by a factor that decays exponentially with temperature. The resulting decrease in the KIE with temperature has been observed both experimentally17 and theoretically17, 27 for C-H activation catalyzed by the soybean lipoxygenase enzyme, which involves a weak hydrogen bond in the low-frequency regime.

Finally, in this analysis, a specific (µ/ν) pair of reactant/product vibronic states is defined to be in the inverted region when ΔGμν0>λ . For PCET reactions, a plot of the rate constant as a function of ΔG0 is not expected to exhibit a turnover as ΔG0 becomes more negative, with a maximum at ΔGμν0=λ , as predicted by Marcus theory for electron transfer.42 Instead, this PCET theory predicts that the rate constant will increase as ΔG0 becomes more negative because excited vibronic product states associated with low free energy barriers and relatively large vibronic couplings become accessible.43 Although a related theory for nonadiabatic proton transfer also predicts the absence of a turnover in the inverted region,20, 21 such turnover behavior has been observed experimentally for nonadiabatic proton transfer.44 A possible explanation for this inconsistency between theory and experiment is the non-linear geometry for proton transfer in the experimental system, requiring a multidimensional treatment of the hydrogen wavefunction.44 These issues are not directly relevant to the present work, which focuses on the temperature dependence of the KIE for a fixed driving force. The dominance of the (0/1) pair of reactant/product vibronic states for this particular driving force, however, is consistent with the theoretical prediction of an absence of turnover behavior in the inverted region. If ΔG0 were to become somewhat more negative, this theory predicts that the (0/2) pair of reactant/product vibronic states would become the dominant contributor. Moreover, this PCET theory predicts that the ln[KIE] will have a maximum near zero driving force and decrease significantly as the driving force becomes more positive or negative because the contributions from excited vibronic states increase as the reaction becomes more asymmetric, and contributions from excited vibronic states decrease the KIE.43 In principle, such predictions could be tested experimentally with chemically modified quinols.

IV. Conclusions

In this paper, we studied the PCET reactions in biomimetic model systems for ubiquinol and plastoquinol oxidation. Our calculations provide a possible explanation for the experimental observation that the apparent activation free energy is greater for hydrogen than for deuterium, leading to an increase in the KIE with increasing temperature, for ubiquinol oxidation. This unusual temperature dependence of the KIE is possible if the following conditions are satisfied: (1) high-frequency proton donor-acceptor motion, which corresponds to a stiff hydrogen bond; (2) small inner-sphere and solvent reorganization energy, so that –ΔG0 > λ for the (0/0) pair of reactant/product vibronic states; (3) the (0/0) pair of reactant/product vibronic states is in the inverted Marcus region, while the (0/1) pair of reactant/product vibronic states is in the normal Marcus region; (4) the (0/1) pair of reactant/product vibronic states is the dominant contributor to the overall rate, which could occur because of the larger overlap between the reactant and product proton vibrational wavefunctions for excited states. If these conditions are satisfied, the free energy barrier for the dominant transition is lower for deuterium than for hydrogen because of the smaller splittings between the vibronic energy levels for deuterium, and the KIE increases with increasing temperature.

Our calculations indicate that the biomimetic UQH2 and PQH2 oxidation reactions in acetonitrile have relatively small reorganization energies and high-frequency proton donor-acceptor modes. For both systems, the (0/0) pair of reactant/product vibronic states is in the inverted Marcus region, while the (0/1) pair of reactant/product vibronic states is in the normal Marcus region. In the UQH2 system, the (0/1) pair of reactant/product vibronic states is the dominant contributor to the overall rate for both hydrogen and deuterium, and the explanation for the unusual inverse temperature dependence of the KIE given above is applicable. In the PQH2 system, however, other pairs of reactant/product vibronic states contribute significantly to the overall reaction rate, leading to more complex temperature dependence of the KIE. As a result, the KIE increases with increasing temperature for UQH2 and decreases with increasing temperature for PQH2. This qualitatively different behavior for the two systems is due to the rather subtle differences in the driving forces and the proton potential energy curves, which impact the free energy barriers and vibronic couplings for each pair of reactant/product vibronic states and thereby determine the relative contributions from each pair.

Although the temperature dependence of the KIE for quinol oxidation is not directly relevant to the function of biological systems such as the cyt bc1 complex, studies of the temperature dependence of the KIE may reveal the underlying mechanism of these reactions. The present model calculations, in conjunction with the experimental observation of the unusual inverse temperature dependence of the KIE for ubiquinol oxidation in the cyt bc1 complex,3 indicate that this reaction occurs in the inverted Marcus region and involves significant contributions from excited vibronic states. These types of fundamental insights may be useful for harnessing and enhancing biological energy transduction processes.

Supplementary Material

1_si_001

Supporting Information Available:

Structures of the truncated gas phase models, truncated solvated models, the full system gas phase models, and the superimposed full system models; partial charges used for the calculation of solvent reorganization energy and proton potential energy curves; molecular mechanical parameters used to fit the gas phase proton potential energy curves; analog of Table 3 at 320 K; estimation of the driving forces; complete reference 29. This material is available free of charge via the Internet at http://pubs.acs.org.

Acknowledgments

We are grateful for many helpful discussions with David Kramer and Jon Cape. We also thank Jon Cape for providing us with his estimates of the driving forces for the two quinol oxidation reactions, as given in Supporting Information. We acknowledge financial support from NSF Grant No. CHE-07-49646 and NIH grant GM56207.

References

1. Cramer WA, Knaff DB. Energy transduction in biological membranes: A textbook of bioenergetics. New York: Springer-Verlag; 1990.
2. Cape JL, Bowman MK, Kramer DM. Trends in Plant Science. 2006;11:46–55. [PubMed]
3. Cape JL, Bowman MK, Kramer DM. J. Am. Chem. Soc. 2005;127:4208–4215. [PubMed]
4. Bell RP. The Tunnel Effect in Chemistry. London: Chapman and Hall; 1980.
5. Pu J, Gao J, Truhlar DG. Chemical Reviews. 2006;106:3140–3169. [PubMed]
6. Nagaoka S-i, Inoue M, Nishioka C, Nishioku Y, Tsunoda S, Ohguchi C, Ohara K, Mukai K, Nagashima U. J. Phys. Chem. B. 2000;104:856–862.
7. Nagaoka S-i, Nishioku Y, Mukai K. Chem. Phys. Lett. 1998;287:70–74.
8. Cukier RI, Nocera DG. Annual Reviews of Physical Chemistry. 1998;49:337–369. [PubMed]
9. Mayer JM. Annu. Rev. Phys. Chem. 2004;55:363–390. [PubMed]
10. Huynh MH, Meyer TJ. Chemical Reviews. 2007;107:5004–5064. [PMC free article] [PubMed]
11. Yamamoto T, Kato S. J. Chem. Phys. 2007;126:224514. [PubMed]
12. Soudackov A, Hammes-Schiffer S. J. Chem. Phys. 1999;111:4672–4687.
13. Soudackov A, Hammes-Schiffer S. J. Chem. Phys. 2000;113:2385–2396.
14. Soudackov A, Hatcher E, Hammes-Schiffer S. J. Chem. Phys. 2005;122:014505. [PubMed]
15. Hammes-Schiffer S, Soudackov AV. J. Phys. Chem. B. 2008;112:14108–14123. [PMC free article] [PubMed]
16. Pu J, Ma S, Gao J, Truhlar DG. J. Phys. Chem. B. 2005;109:8551–8556. [PubMed]
17. Knapp MJ, Rickert KW, Klinman JP. J. Am. Chem. Soc. 2002;124:3865–3874. [PubMed]
18. Meyer MP, Klinman JP. Chem. Phys. 2005;319:283–296. [PMC free article] [PubMed]
19. Siebrand W, Smedarchina Z. J. Phys. Chem. B. 2004;108:4185–4195.
20. Kiefer PM, Hynes JT. J. Phys. Chem. A. 2004;108:11793–11808.
21. Kiefer PM, Hynes JT. J. Phys. Chem. A. 2004;108:11809–11818.
22. Antoniou D, Schwartz SD. J. Phys. Chem. B. 2001;105:5553–5558.
23. Hay S, Pudney C, Hothi P, Johannissen LO, Masgrau L, Pang J, Leys D, Sutcliffe MJ, Scrutton NS. Biochem. Soc. Trans. 2008;36:16–21. [PubMed]
24. Marcus RA, Sutin N. Biochim. Biophys. Acta. 1985;811:265–322.
25. Georgievskii Y, Stuchebrukhov AA. J. Chem. Phys. 2000;113:10438–10450.
26. Skone JH, Soudackov AV, Hammes-Schiffer S. J. Am. Chem. Soc. 2006;128:16655. [PubMed]
27. Hatcher E, Soudackov AV, Hammes-Schiffer S. J. Am. Chem. Soc. 2007;129:187–196. [PubMed]
28. Cape J. personal communication. 2009. In.
29. Frisch MJ, et al. Gaussian 03, revision C.03. Pittsburgh, PA: Gaussian, Inc.; 2003.
30. Lee C, Yang W, Parr PG. Physical Review B. 1988;37:785–789. [PubMed]
31. Becke AD. J. Chem. Phys. 1993;98:5648–5652.
32. Stephens PJ, Devlin FJ, Chablowski CF, Frisch MJ. J. Phys. Chem. 1994;98:11623.
33. Hay PJ, Wadt WR. J. Chem. Phys. 1984;82:270.
34. Miertus S, Scrocco E, Tomasi J. Chem. Phys. 1981;55:117–129.
35. Barone V, Cossi M, Tomasi J. Journal of Chemical Physics. 1997;107:3210.
36. Breneman CM, Wiberg KB. J. Comput. Chem. 1990;11:361–373.
37. Basilevsky MV, Rostov IV, Newton MD. Chem. Phys. 1998;232:189–199.
38. Newton MD, Basilevsky MV, Rostov IV. Chem. Phys. 1998;232:201–210.
39. Iordanova N, Decornez H, Hammes-Schiffer S. J. Am. Chem. Soc. 2001;123:3723–3733. [PubMed]
40. Biner M, Burgi H-B, Ludi A, Rohr C. J. Am. Chem. Soc. 1992;114:5197–5203.
41. Ishikita H, Soudackov AV, Hammes-Schiffer S. J. Am. Chem. Soc. 2007;129:11146–11152. [PubMed]
42. Marcus RA. J. Chem. Phys. 1989;93:3078–3086.
43. Edwards SJ, Soudackov AV, Hammes-Schiffer S. Journal of Physical Chemistry A. in press
44. Peters KS. Acc. Chem. Res. 2009;42:89–96. [PubMed]