Adiabaticity of the Proton-Coupled Electron-Transfer Step in the

Apr 7, 2015 - Structure of the Active-Site Structure of NiSOD and the O2– .... All electronic structure calculations were performed using the ORCA s...
0 downloads 0 Views 1MB Size
Download PDF
Article pubs.acs.org/JPCB

Adiabaticity of the Proton-Coupled Electron-Transfer Step in the Reduction of Superoxide Effected by Nickel-Containing Superoxide Dismutase Metallopeptide-Based Mimics Jason Shearer,* Jennifer C. Schmitt, and Heather S. Clewett Department of Chemistry, University of Nevada, 1664 North Virginia Street, Reno, Nevada 89557, United States S Supporting Information *

ABSTRACT: Nickel-containing superoxide dismutases (NiSODs) are bacterial metalloenzymes that catalyze the disproportionation of O2−. These enzymes take advantage of a redox-active nickel cofactor, which cycles between the Ni(II) and Ni(III) oxidation states, to catalytically disprotorptionate O2−. The Ni(II) center is ligated in a square planar N2S2 coordination environment, which, upon oxidation to Ni(III), becomes five-coordinate following the ligation of an axial imidazole ligand. Previous studies have suggested that metallopeptide-based mimics of NiSOD reduce O2− through a proton-coupled electron transfer (PCET) reaction with the electron derived from a reduced Ni(II) center and the proton from a protonated, coordinated NiII−S(H+)−Cys moiety. The current work focuses on the O2− reduction half-reaction of the catalytic cycle. In this study we calculate the vibronic coupling between the reactant and product diabatic surfaces using a semiclassical formalism to determine if the PCET reaction is proceeding through an adiabatic or nonadiabatic proton tunneling process. These results were then used to calculate H/D kinetic isotope effects for the PCET process. We find that as the axial imidazole ligand becomes more strongly associated with the Ni(II) center during the PCET reaction, the reaction becomes more nonadiabatic. This is reflected in the calculated H/D KIEs, which moderately increase as the reaction becomes more nonadiabatic. Furthermore, the results suggest that as the axial ligand becomes less Lewis basic the observed reaction rate constants for O2− reduction should become faster because the reaction becomes more adiabatic. These conclusions are in-line with experimental observations. The results thus indicate that variations in the axial donor’s ability to coordinate to the nickel center of NiSOD metallopeptide-based mimics will strongly influence the fundamental nature of the O2− reduction process.

■

INTRODUCTION Nickel-containing superoxide dismutase (NiSOD) is the most recently discovered SOD and is widely expressed in soil and aquatic bacteria.1−7 As is the case with all SODs, NiSOD catalyzes the disproportionation of superoxide (O2−) into hydrogen peroxide and dioxygen by cycling between reduced and oxidized oxidation states, NiII and NiIII in the case of NiSOD Ni IISOD + 2H+ + O2− → Ni IIISOD + H 2O2

(1)

Ni IIISOD + O2− → Ni IISOD + O2

(2)

Scheme 1. Structure of the Active-Site Structure of NiSOD and the O2− Disproportionation Reaction Catalyzed by NiSOD

used this convenience of nature to construct NiSOD metallopeptide-based mimics (NiSOD MBMs) derived from the NiSOD primary sequence.1,11−15 Many of these NiSOD MBMs are catalytically active, and in some cases can display catalytic rates approaching that of the metalloenzyme.1 A number of experimental studies have been performed aimed at understanding the reactivity and mechanism of these NiSOD MBMs. One study demonstrated that the axial ligand

The nickel center in reduced NiSOD is contained within a square planar NiN2S2 coordination environment with ligands derived from the N-terminal amine nitrogen atom, the amidate nitrogen atom from Cys(2), and sulfur atoms from Cys(2) and Cys(6). Oxidation to the NiIII oxidation state yields a square pyramidal NiN3S2 nickel center following the ligation of the His(1) imidazole δ-nitrogen atom (Scheme 1).8−10 Thus, all of the ligating residues to nickel are found within the first six residues from the NiSOD N-terminus. We, and others, have © 2015 American Chemical Society

Received: March 19, 2015 Published: April 7, 2015 5453

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B

Chart 1. General active-site structure of reduced NiSOD MBMs (A) and the structure of computational models 1 (B), 2 (C), and 2-off (D).

Scheme 2. Proposed Mechanism for the Reduction of Superoxide Facilitated by NiSOD MBMs

subsequently applied to PCET reactions by Hammes-Schiffer and coworkers.22 In addition to probing the adiabaticity of the superoxide reduction reaction facilitated by these NiSOD MBMs, we will explore the role, if any, that the Lewis basicity of the axial imidazole ligand to Ni(II) center has on the nature of the PCET process.

to nickel likely remains associated with the Ni(II) center under rapid catalytic conditions.12 It is only under slow turnover conditions (i.e., low superoxide loads) that the axial imidazole ligand dissociates from the Ni(II) center. Another study demonstrated that the nature of the imidazole ligand had an influence on the rate of O2− disproportionation catalysis; when the axial ligand was made less Lewis basic, the rate of catalysis increased by two orders of magnitude.15 In two more recent studies the mechanism by which these NiSOD MBMs disproportionate O2− was probed.16,17 In one study examining the superoxide disproportionation kinetics effected by a NiSOD MBM with an unnatural εN-methylhistidine (Ni-m1Me, Chart 1), it was speculated that superoxide reduction likely proceeds through a proton-coupled electrontransfer (PCET) reaction.17 This was evidenced, in part, by a large solvent H/D kinetic isotope effect (KIE; room temperature KIE = 20(4)). It was also proposed that the proton and electron were derived from a protonated-coordinated NiII− S(H+)−Cys moiety. A follow-up study pointed toward the Cys(6) sulfur atom as the likely candidate for the protonated cysteinate in the NiII−S(H+)−Cys moiety.16 Furthermore, it was suggested that one water molecule was contained within a hydrogen bonding network to the Cys(2) and Cys(6) sulfur atoms, which thermodynamically aided the overall O 2 − reduction process (Chart 1, Scheme 2). Lastly, computational studies suggested that the PCET reaction proceeded via a hydrogen atom transfer reaction (HAT reaction, i.e., an electronically adiabatic proton tunneling process). Although this mechanism is not relevant to the metalloenzyme itself,18 it is an interesting PCET process that may be of fundamental importance in the PCET mechanisms of other redox-active thiolate-ligated transition-metal systems. This study aims at evaluating the PCET reaction facilitated by NiSOD MBMs in a more rigorous manner than that previously employed. The conclusion that O2− reduction proceeded through a HAT process was based solely on MO arguments. Conclusions regarding the mechanism of PCET reactions can be difficult to make on the basis of such arguments alone.19,20 Therefore, we further probed the adiadibicity of the reaction using a semiclassical analysis derived by Georgievskii and Stuchebrukhov, 21 which has been

■

THEORETICAL MODEL In this study, we examine the transfer of a proton and an electron from a Ni(II)−S(H+)−R moiety to a superoxide anion Ni(II)−S(H+)−R + O2− → Ni(III)−S−R + HO2−

(3)

The following treatment utilizes a two-state valence bond model to approximate the reactant and product diabatic surfaces for the net hydrogen transfer reaction. As previously demonstrated by Hammes-Schiffer,22 the energies of these two surfaces, represented as Hda,1(q) and Hda,2(q), are a function of proton position, q. The electronic coupling between the surfaces, VET, is assumed to be a constant with relation to the proton position at a given proton donor−acceptor distance but will vary from one donor−acceptor distance to another. Hda,1(q), Hda,2(q), and VET are used to generate the following valence bond Hamiltonian matrix ⎡ H (q) V ET ⎤ ⎢ da,1 ⎥ ⎢ ET ⎥ Hda,2(q)⎦ ⎣ V

(4)

Diagonalization of the valence bond Hamiltonian matrix yields the electronic adiabatic surfaces U1(q) =

U2(q) =

5454

1 (Hda,1(q) + Hda,2(q)) 2 1 − (Hda,1(q) − Hda,2(q))2 + 4(V ET)2 2

(5)

1 (Hda,1(q) + Hda,2(q)) 2 1 + (Hda,1(q) − Hda,2(q))2 + 4(V ET)2 2

(6)

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B where U1(q) and U2(q) are the ground- and excited-state adiabatic surfaces, respectively. The calculation of the proton-tunneling matrix elements, Tx, follows that outlined by Georgievskii and Stuchebrukhov21 and Hammes-Schiffer.22 In brief, in the nonadiabatic limit, the tunneling matrix element between the reactant and product diabatic surfaces, Tna, is related to the Franck−Condon overlap of the reactant and product proton vibrational wave functions, φR and φP, by T na = V ET⟨φR |φ P⟩

and reactant vibronic wave functions. Thus, for a given donor− acceptor distance displaced from the equilibrium donor− acceptor distance, Tx will vary by T x = Tox exp[−α(R − R o)]

where Txo is the tunneling matrix element at the equilibrium geometry.23,24 The distance dependence of tunneling matrix elements can be used to approximate KIEs. Previous work has shown that for nonadiabatic transitions between the ground -tate vibronic wave functions the H/D kinetic isotope effect (KIE) can be approximated using

(7)

which results in a relatively small tunneling element. In the adiabatic limit, the resulting tunneling matrix element Tad is larger than Tna, although it is still considerably smaller than VET. In this manuscript, we use the following quasiclassical approximation for the calculation of Tad T ad

⎡ 1 ≈ 0.17ℏ νDνA exp⎢ − ⎣ ℏ

∫q

qB

A

KIE ≈

(9)

where the prefactor κ is defined as exp[p ln p − p] Γ(p + 1)

(10)

In eq 10, Γ(x) is the Gamma function, and the adiabticity parameter, p, is defined as the ratio of the proton tunneling time (τp) over the electronic transition time (τe) τp p= τe (11) The values τp and τe are related to VET by

τp ≈

V ET |ΔF |νt

(12)

ℏ V ET

(13)

and τe ≈

where |ΔF| is the difference in the slopes at the crossing point of the reactant and product surfaces and νt is the proton tunneling velocity νt =

2(Ec − Et) m

⎡ 2k T ⎤ exp⎢ B 2 (αH2 − αD2)⎥ ⎣ MΩ ⎦

(16)

COMPUTATIONAL METHODS All electronic structure calculations were performed using the ORCA software package (versions 2.9 and 3.0.2).26 Three different computational models of the NiSOD MBMs were considered, one containing a glycine moiety (1) and the others containing an εN-methyl histidine moiety with the imidazole δnitrogen atom either associated with (2) or not coordinated to the Ni(II)-center (2-off, Chart 1). The initial input geometry for 1 was taken from the reactant state of the O2− + 1 →HO2− + 1ox reaction previously investigated. The initial input geometry for 2 was taken from 1 with the glycine H group replaced with a −CH2−εN-methyl imidazole moiety. Following initial geometry optimizations at the BP86/def2-tzvp level of theory, final geometry optimizations were performed at the B2PLYP/def2-tzvp level of theory using the COSMO solvation model (ε = 80.4; r = 1.33).27−29 All electronic structure calculations were performed using a doublet ground state. To obtain transition-state structures, we constructed an initial transit of the hydrogen atom from the thiolate sulfur (donor) to the superoxide oxygen (acceptor) atoms in 1, 2, and 2-off by moving the transferred proton on an evenly spaced grid of 20 hydrogen atom positions between the donor and acceptor atoms. The subsequent transition-state structures were then optimized starting from a geometry at the extrapolated maximum from the relaxed PES scan. The vibrational analyses of the transition-state structures indicated one imaginary frequency corresponding to the S−H−O vibrational mode (ν =−502i cm−1 for 1, −584i cm−1 for 2, and −541i cm−1 for 2off). The input quantities used to determine the tunneling matrix elements for the 1 + O2− and 2(-off) + O2− reactions utilized the previously described transition-state structures. All atoms, except for the transferring H atom, and the donor (S) and acceptor (O) distances were fixed. The transferring H atom was then moved along a systematic grid between the donor and acceptor atoms. Single-point energies of different transferring hydrogen atom positions were determined along a 3D grid of 75 points centered about the vector describing the S−H−O vibrational mode at the transition-state equilibrium geometries. The minimum energy pathway for the transferring hydrogen atom along this 3D grid was then determined and used to construct the transit vector for the transfer of the hydrogen atom from the donor to acceptor atoms. This transit vector

where νD/A are the vibrational frequencies in the donor and acceptor wells, qA/B are the turning points of the classical trajectory at the energy E, m is the mass of the proton, and U(q) is the potential energy as a function of proton position, q.21,25 The semiclassical expression for the tunneling matrix element, Tsc, for the whole range of variations in VET is expressed as

2πp

x 2 |To,D |

■

(8)

κ=

x 2 |To,H |

where M and Ω are the effective mass and frequency associated with the R mode, T in the exponential term is temperature, and kB is the Boltzmann constant.24

⎤ 2m[U (q) − E] dq⎥ ⎦

T sc = κT ad

(15)

(14)

with Ec corresponding to the crossing point energy, Et corresponding to the tunneling energy, and m corresponding to the mass of the proton. To a first approximation, the dependence of the tunneling matrix element on donor−acceptor distance, R, should be exponential because of the distance dependence of the product 5455

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B corresponded to the vector describing the S−H−O imaginary vibrational mode. Single-point calculations for the H-atom transfer reaction at the O2− reduction reaction transition-state equilibrium geometry were then performed for different hydrogen atom positions, q, along this transit vector generating a 1D PES at the NEVPT2 level of theory.30 The NEVPT2 calculations utilized a CAS-SCF(3,6) active space and the tzvp basis set. Natural orbitals from MP2/tzvp level calculations were used as input orbitals for the NEVPT2 calculations. The resulting 1D PESs were constructed from 35 individual single point calculations with a tighter Δq grid in the region of the transition state. To calculate the dependence of the tunneling matrix elements, T, on donor/acceptor distance, we generated 1D PESs with different donor−acceptor distances in an analogous manner as that previously outlined. These calculations were performed at the B2PLYP/tzvp level of theory because this level yields results similar to the NEVPT2 level calculations for these systems. The values used for Ω and M in eq 16 were determined as previously outlined.32 The calculated values used for VET were taken as half the energy separating the ground- and excited-state surfaces at the TS hydrogen atom position on each corresponding PES. The reactant and product diabatic surfaces displayed in eqs 5 and 6 were approximated as Morse potentials with an exponential wall function at large donor/acceptor distances. Vibrational eigenvalues/eigenfunctions were calculated using the Fourier Grid Hamiltonian method.31

■

Figure 2. Ground (blue) and excited (red) state adiabatic PESs for the transfer of a hydrogen atom from 1 to O2− (A) and 2 to O2− (B). The distance r = 0 Å indicates the hydrogen atom position at the transition state with negative distances closer to the donor S atom and positive distances to the acceptor O atom. The red and blue circles represent the calculated energies along the adiabatic PESs. The dashed green lines represent the reactant and product diabatic surfaces. For these PESs, the donor/acceptor distance for 1···O2− is 3.070 Å and 2···O2− is 3.027 Å.

RESULTS Adiabadicity and Tunneling Matrix Elements for the O2− Reduction Reaction Facilitated by 1 at the Equilibrium Transition-State Geometry. The transitionstate structure for the 1 + O2− → 1ox + HO2− reaction is shown in Figure 1. As previously determined for 1, the structure of the

primarily a S(3p)−Ni(3dπ)* orbital with virtually no AO character derived from the protonated sulfur atom. In addition, there is minimal contributions from the O2− π* orbital fragment that is contained in the molecular xy plane. Furthermore, there is also some water O 2px contribution to the SOMO. As the hydrogen atom progresses to its transitionstate position, the amount of S(3p)−Ni(3dπ) character decreases while O2− π* and water O 2px character increase. Also, the SOMO at the transition state now contains a significant contribution from the protonated sulfur atom’s S(3p) orbital that is antibonding to the O2− π* MO fragment to the SOMO. At the product state, the SOMO is now predominantly composed of a “protonated” S(3p)−O2− π* antibonding interaction.33 Figure 2a displays the NEVPT2 ground- and excited-state adiabatic surfaces for the net hydrogen atom transfer at the transition-state geometry along with reactant and product diabatic surfaces. The electronic coupling, VET, was taken as half the energy difference between the ground- and excited-state adiabatic surfaces at the TS. We find reasonably strong coupling between the surfaces with a VET = 8360 cm−1. Application of eqs 7 and 8 yields tunneling matrix elements at the nonadiabatic (Tna) and adiabatic limits (Tad) of 6.41 and 34.1 cm−1, respectively.

Figure 1. From left to right: Transition-state geometry of 1 and the highest occupied molecular orbital at the reactant, product-state geometries on the diabatic surfaces, and the TS geometry along the adiabatic surface.

1···O2− complex at the transition state contains the O2− molecule within the amine H−N−Ni(II)−S−H cleft with the transferred hydrogen atom between the sulfur atom and O2− oxygen atom. In this structure, the proton donor atom (S) is situated 3.070 Å from the proton acceptor (O) atom. An initial 1D potential energy surface for the reduction of O2− by 1 was generated at the NEVPT/CAS-SCF(3,6) level of theory by moving the transferred hydrogen atom between the donor and acceptor positions and keeping all of the other atoms fixed at their transition state equilibrium positions. The resulting PES is displayed in Figure 2. We note that because 1 possesses a low degree of multiconfigurational character similar results are obtained at the B2PLYP level of theory.34 Figure 1 displays iso-surface plots of the SOMO of 1···O2− at the transition-state geometry and minima on the reactant and product diabatic surfaces. At the reactant state, the SOMO is 5456

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B Hammes-Schiffer and coworkers have used the ratio of the proton tunneling time (τp) over the electronic transition time (τe) as an indication of the overall adiabaticity of the protontunneling reaction; this ratio is the adiabaticity parameter p.19,20,22 As the proton tunneling time becomes much longer than the electron-transition time (i.e., when p ≫ 1), the electrons can relax in response to the movement of the proton, and the proton-tunneling reaction is described as an adiabatic process. In contrast, when the proton-tunneling time is shorter than the electronic-transition time (i.e., when p ≪ 1), the electrons do not have sufficient time to relax, and the protontunneling reaction is described as a nonadiabatic process. For the 1 + O2− PCET reaction the calculated values for τe (0.668 fs) and τp (0.122 fs) result in a value for p = 0.182. The corresponding value for p is neither exceedingly large nor small, indicating the proton-tunneling reaction could be described as adiabatic with significant nonadiabatic character. The resulting values for κ and Tsc support this supposition. Application of the derived value for p to eq 10 yields κ = 0.709 and a resulting value for Tsc = 24.2 cm−1. Adibadicity and Tunneling Matrix Elements for the O2− Reduction Reaction Facilitated by 2 at the Equilibrium Transition-State Geometries. Two transitionstate structures for the reduction of O2− by the histidinecontaining computational model were explored. In one case, 2off, the N-methyl histidine imidazole ligand is positioned over the Ni(II) center, with the imidazole ring perpendicular to the molecular z axis. In the second case, 2, the imidazole δ-nitrogen is associated with the Ni(II) center. As was previously demonstrated, there is a shallow local minima corresponding to the associated Ni(II)−Nimidazole structure (2), which lies 18.1 kJ/mol higher in energy than 2-off at the B2PLYP/def2-tzvp level of theory.12 In the case of 2, the axial Ni−N distance is considerably longer than what would be expected for a Ni−N bond (2.623 Å vs ∼2.0 Å) and is best considered as a Ni−N interaction as opposed to a Ni−N bond. At the 2···O2− transition-state equilibrium geometry, there is a contraction of the proton donor (S) acceptor (O) distance relative to that obtained for 1···O2−, with an equilibrium S−O bond length of 3.027 Å. The shorter equilibrium donor−acceptor bond length is likely the result of the higher degree of negative charge on the superoxide anion at the transition state of 2···O2− versus 1··· O2−. Not surprisingly, the results obtained for 2-off···O2− were similar to the results obtained for 1···O2−. This might have been expected considering the lack of axial interaction between the His-imidazole ligand and the Ni(II) center. Because of the similar results obtained between the two models, 2-off···O2− was not fully analyzed and will not be further discussed. In contrast, the results obtained for 2···O2− were considerably different than those obtained for 1···O2−, owing to the influence of the histidine imidazole on the electronic structure of the Ni center. Figure 3 displays the SOMO at the transition-state geometry and minima on the reactant and product diabatic surfaces for the superoxide reduction reaction facilitated by 2···O2− at the NEVPT2/CAS-SCF(3,6) level of theory. The reactant-state wave function is similar to that observed for 1···O2−, which is primarily a S(3p)−Ni(3dπ)* antibonding MO with some contribution from the amidate N(2p) orbital and water O(2p) orbitals. The primary differences between 1···O2− and 2···O2− are (a) there is now significantly more protonated sulfur S(3p) character to the wave function and (b) there is an additional

Figure 3. From left to right: Transition-state geometry of 2 and the highest occupied molecular orbital at the reactant, product-state geometries on the diabatic surfaces, and the TS geometry along the adiabatic surface.

contribution from the axial imidazole δN(2p) orbital to the wave function. At the transition-state hydrogen atom position, the SOMO is qualitatively similar to the reactant state SOMO with slightly more superoxide O(2p) and protonated S(3p) contributions to the wave function at the expense of S(3p)− Ni(3dπ)* antibonding character. As was found with 1···O2−, at the product state, the SOMO for 2···O2− contains significantly more O(2p) character derived from the O2 moiety than found at the transition state; however, unlike the case of 1···O2−, where the wave function is almost exclusively O2−S(3p)* in character, the product-state wave function for 2···O2− displays a significant contribution from the NiN3S2 moiety. This seems to suggest that the PCET O2− reduction process facilitated by 1 is fundamentally different than that facilitated by 2. This is supported by the analysis of the PCET PES for the superoxide reduction process in 2···O2−. Figure 2b displays the NEVPT2/CAS-SCF(3,6) generated ground- and excited-state surfaces for the 2 facilitated O2− reduction reaction. The separation between the ground- and excited-state surfaces at the transition-state hydrogen atom position is considerably smaller than that obtained for 1···O2−, translating into a VET = 2100 cm−1. The lower degree of electronic coupling between the two diabatic surfaces results in a considerably longer electronic transition time than that obtained for 1···O2− (τe= 2.62 vs 0.668 fs). In contrast, the proton-tunneling time remains comparable to that obtained in the analysis of 1···O2− (τp= 0.155 vs 0.122 fs). Application of eqs 10 and 11 yields p = 0.0592 and κ = 0.502. This results in a value for Tsc= 10.8 cm−1, which is closer to the value obtained for the nonadiabatic tunneling matrix element (Tna = 7.22 cm−1) than the corresponding adiabatic tunneling matrix element (Tad = 21.6 cm−1). This suggests that proton tunneling in the 2 + O2− reaction is primarily a nonadiabatic process with some adiabatic character. Calculated H/D Kinetic Isotope Effects For O 2 − Reduction by 1 and 2. We sought to gain insight into the nature of the experimentally observed H/D KIE values obtained for superoxide disproportionation by these NiSOD MBMs. This process involves determining the influence of donor/acceptor distance on Tx. PESs were generated as previously described using a range of donor/acceptor distances, R, with all other atoms fixed at their equilibrium positions. The resulting distance dependence of the tunneling matrix elements was then used to estimate the H/D KIE for these reactions according to eqs 15 and 16 (Tables 1 and 2). These results are based on PESs generated from B2PLYP/tzvp calculations because such results compare well with the more costly NEVPT2 calculations. The distance dependence of Tx for the PCET reaction facilitated by 1 and 2 is displayed in Figure 4, the corresponding values for αH/D are given in Tables 1 and 2, and the calculated 5457

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B

Table 1. Dependence of Tx (in cm−1) on Donor (S) Acceptor (O) Distance (R) and the Corresponding Values for α for the Superoxide Reduction Reaction Facilitated by 1a R = 2.959 Å R = 3.019 Å Ro = 3.070 Å R = 3.119 Å R = 3.169 Å α a

Tad H

TscH

Tna H

Tad D

TscD

Tna D

272 98.6 34.1 20.8 4.94 18

172 61.6 24.2 12.9 3.04 18

153 14.1 6.41 2.39 0.82 19

58.5 22.1 4.72 2.34 1.05 21

46.2 14.0 2.95 1.46 0.64 22

12.2 5.11 0.87 0.11 0.01 31

Values at Ro were determined at the NEVPT2/CAS-SCF(3,6) level of theory. All other values were determined at the B2PLYP level of theory.

Table 2. Dependence of Tx (in cm−1) on Donor (S) Acceptor (O) Distance (R) and the Corresponding Values for α for the Superoxide Reduction Reaction Facilitated by 2a R = 2.988 Å R = 2.995 Å Ro = 3.027 Å R = 3.051 Å R = 3.063 Å α

Tad H

TscH

Tna H

Tad D

TscD

Tna D

52.3 45.3 21.6 14.2 10.6 21

29.5 25.9 10.8 7.14 5.42 23

21.9 13.2 7.22 3.86 2.71 25

8.28 7.60 2.69 1.93 1.56 24

5.64 4.31 1.63 1.09 0.89 26

4.89 3.91 1.17 0.74 0.53 30

Table 3. Calculated H/D KIE Values for the Superoxide Reduction Reaction Facilitated by 1 and 2 Using the SemiClassical Tunneling Matrix Element and Tunneling Matrix Elements at the Adiabatic and Non-Adiabatic Limits 1 2

Tsc

Tad

Tna

22 30

27 46

1.7 19

KIE values over the temperature range of 273−373 K regardless of the tunneling matrix element employed, as is observed experimentally.

a

Values at Ro were determined at the NEVPT2/CAS-SCF(3,6) level of theory. All other values were determined at the B2PLYP level of theory.

■

DISCUSSION Utilizing a computational/theoretical treatment, we have probed the reduction of superoxide by NiSOD MBMs. On the basis of two recent studies this reaction is postulated to occur via a PCET process. Two different computational models were considered. One of these models contains an N-terminal alanine residue, and the other model contains an N-terminal ε N-methyl-histidine residue. In both cases, a PCET mechanism for O2− reduction appears to be a viable mechanism. For both four-coordinate 1 and quasi-five-coordinate 2, the reduction of

H/D KIE values are given in Table 3. In the case of both 1··· O2− and 2···O2−, the calculated H/D KIE values are largest when Tad is utilized, followed by Tsc and then Tna. Furthermore, the calculated H/D KIEs are all moderately larger in the case of 2···O2− than 1···O2−. For 1···O2−, the resulting KIEs range from 27 (Tad) to 1.7 (Tna), with the semiclassical term yielding a KIE of 22. These increase to 46 (Tad), 19 (Tna), and 30 (Tsc) for the 2···O2− facilitated reaction. We note that for these models there is virtually no temperature dependence on the calculated H/D

Figure 4. Dependence of Tx on donor (S) acceptor (O) distance for the hydrogen and deuterium transfer reaction facilitated by 1 (left) and 2 (right). The circles and solid trend lines depict the results for the hydrogen transfer reaction, while the triangles and dashed trend lines depict results for the deuterium transfer reaction. The top panels display the results for Tad, the middle panels display the results for Tsc, and the bottom panels display the results for Tna. 5458

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B O2− by the Ni(II)−S(H+)−Cys moiety appears to proceed by thermodynamically favorable processes with no intermediates formed during the reaction. The two models probed can be seen as extremes in terms of axial donation to the Ni(II) center, ranging from no axial donation to a strongly Lewis basic imidazole ligand. This is an important aspect of O2− reduction by these NiSOD MBMs. In one study examining the influence of axial imidazole Lewis basicity on the properties of NiSOD MBMs, it was demonstrated that as the axial imidazole ligand became more electron rich, the rate of the catalytic superoxide disproportion reaction decreased by nearly two orders of magnitude.15 The NiSOD MBM containing the least electron-rich imidazole donor is in-fact so weakly donating that it would be adequately approximated by 4-coordinate 1. In contrast, the NiSOD MBM with the most electron-rich axial ligand contained the εNmethyl imidazole and is thus approximated by 2. In another study probing the structural rearrangement taking place about the Ni-center during the Ni(II)/Ni(III) redox process, it was discovered that a quasi-five coordinate Ni(II) species would be the dominant reduced form under rapid ET conditions (Scheme 3).12 Such conditions would be encountered during

state μ, λμν is the solvent reorganization energy for states μ and ν, and ΔGμν° is the thermodynamic driving force of the PCET reaction. Owing to the dependence of the rate constant on the square of the tunneling matrix element, we would predict that the weaker the axial interaction, which yields a more adiabatic proton tunneling process (i.e., larger values for the tunneling matrix element), the faster the O2− reduction kinetics. This is what we observed experimentally; the NiSOD MBM with less Lewis basic imidazole ligands yields faster catalytic rates. Other factors, such as differences in the thermodynamic driving force and reorganization energies, will further attenuate the unimolecular rate constant for the PCET reaction. This will be the subject of a future study. We also find a moderate influence on H/D KIEs as the PCET reaction switches from an adiabatic to a nonadiabatic process. For 2, which performs O2− reduction through a predominantly nonadiabatic process, we find that the calculated H/D KIEs increase relative to 1, which affects a predominantly adiabatic PCET O2− reduction reaction. Of the NiSOD MBMs that have been explored, Nim1Me has the most basic imidazole ligand. We would therefore predict that the other catalytically active NiSOD MBMs would all yield smaller solvent H/D KIEs than Nim1Me. We end by making a few general observations concerning the approximations previously used and the general applicability of these methods. The calculation of Tad, and consequently Tsc, is based on a quasi-classical approximation as opposed to more rigorous treatments that have been developed.36 This will undoubtedly lead to errors associated with the calculated rate constants and H/D KIEs when utilized as input quantities. Furthermore, the analysis used to calculate the PESs and subsequent tunneling matrix elements only considered the dominant vibrational mode (the R mode) and ignored the contributions that the other vibrational modes would have on Tx. This approximation may not be entirely valid because it has been shown that consideration of additional vibrational modes in the calculation of Tx can dramatically improve the computational accuracy of Tx when compared with experiment.37−39 Nonetheless, use of these approximations yielded calculated KIE values and predicted rate constant trends that are consistent with experimental results. This seems to suggest that these approximations are qualitatively valid, for this system at least. An analysis of other systems will need to be performed to determine if the previously described results are the product of serendipitous error cancellations or if this treatment is generally applicable to other systems. Another observation concerns the level of theory employed in this study. It might have been predicted that an approach that rigorously considers the multiconfigurational character of the system would be required to accurately capture the electronic structure for a reaction involving the movement of a proton and electron from one molecule to another; however, the superoxide complexes of both 1 and 2 possessed surprisingly little multiconfigurational character. (For both systems, the SOMO contributed ∼90% to the overall wave function.) This is despite the fact that these systems contained a transition-metal ion and superoxide radial. Owing to the low degree of multiconfigurational character of 1···O2− and 2···O2−, PESs constructed at a lower level of theory (the B2PLYP double-hybrid functional) compared well with the CAS-SCF generated PESs. Although this lower level of theory may not be appropriate for other systems, the previously described results suggest that this level of theory yields generally applicable

Scheme 3. Square Scheme Depicting Variable Imidazole Ligation upon Changes in Oxidation State for the Catalytic Reaction Performed by NiSOD MBMs

rapid O2− disproportionation. In fact, on the basis of previously determined rate and equilibrium constants we predict that the four-coordinate species will play an insignificant role in O2− disproportionation catalysis by NiSOD MBMs containing electron-rich axial donors.35 Our findings suggest that as the axial imidazole becomes more strongly interacting with the Ni(II) center, the overall proton tunneling event in the PCET reaction changes from an adiabatic to a nonadiabatic process. This will influence the rate of the overall PCET process and observed KIE values. For a PCET reaction at a fixed donor−acceptor distance, the unimolecular rate constant can be expressed as k=

∑ Pμ ∑ μ

ν

x 2 |Tμν |

ℏ

π × λμνkBT

⎡ (ΔG ° + λ )2 ⎤ μν μν ⎥ exp⎢ − 4λμνkBT ⎢⎣ ⎥⎦

(17)

where the subscripts μ and ν correspond to reactant and product vibrational states, Pμ is the Boltzmann probability for 5459

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B

(10) Herbst, R. W.; Guce, A.; Bryngelson, P. A.; Higgins, K. A.; Ryan, K. C.; Cabelli, D. E.; Garman, S. C.; Maroney, M. J. Biochemistry 2009, 48, 3354. (11) Shearer, J.; Long, L. M. Inorg. Chem. 2006, 45, 2358. (12) Neupane, K. P.; Gearty, K.; Francis, A.; Shearer, J. J. Am. Chem. Soc. 2007, 129, 14605. (13) Tietze, D.; Voigt, S.; Mollenhauer, D.; Tischler, M.; Imhof, D.; Gutmann, T.; Gonzalez, L.; Ohlenschlager, O.; Breitzke, H.; Gorlach, M.; Buntkowsky, G. Angew. Chem., Int. Ed. 2011, 50, 2946. (14) Schmidt, M.; Zahn, S.; Carella, M.; Ohlenschlager, O.; Gorlach, M.; Kothe, E.; Weston, J. ChemBioChem. 2008, 9, 2135. (15) Shearer, J.; Neupane, K. P.; Callan, P. E. Inorg. Chem. 2009, 48, 10560. (16) Shearer, J.; Peck, K. L.; Schmitt, J. C.; Neupane, K. P. J. Am. Chem. Soc. 2014, 136, 16009. (17) Shearer, J. Angew. Chem., Int. Ed. 2013, 52, 2569. (18) Ryan, K. C.; Guce, A. I.; Johnson, O. E.; Brunold, T. C.; Cabelli, D. E.; Garman, S. C.; Maroney, M. J. Biochemistry 2015, 54, 1027. (19) Hammes-Schiffer, S.; Stuchebrukhov, A. Chem. Rev. 2010, 110, 6939. (20) Layfield, J. P.; Hammes-Schiffer, S. Chem. Rev. 2014, 114, 3466. (21) Georgievskii, Y.; Stuchebrukhov, A. A. J. Chem. Phys. 2000, 113, 10438. (22) Skone, J. H.; Soudackov, A. V.; Hammes-Schiffer, S. J. Am. Chem. Soc. 2006, 128, 16655. (23) Hatcher, E.; Soudackov, A.; Hammes-Schiffer, S. Chem. Phys. 2005, 319, 93. (24) Edwards, S. J.; Soudackov, A.; Hammes-Schiffer, S. J. Phys. Chem. A 2009, 113, 2117;Hatcher, E.; Soudackov, A. V.; HammesSchiffer, S. J. Am. Chem. Soc. 2007, 129, 187. (25) Fain, B. Theory of Rate Processes in Condensed Media; Springer: New York, 1980. (26) Nesse, F. The ORCA program system. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 73. (27) (a) Grimme, S. J. Chem. Phys. 2006, 124, 034108;(b) Schwabe, T.; Grimme, S. Phys. Chem. Chem. Phys. 2006, 8, 4398. (28) Ahlrichs, R. Basis sets were obtained from the ftp server of the Karlsruhe quantum chemistry group: ftp.chemie.uni-karlsruhe.de/pub/ bases. (29) Klamt, A.; Schüümann, G. J. Chem. Soc., Perkin Trans. 1993, 2, 799. (30) (a) Angeli, C.; Cimiraglia, R.; Evangelisti, S.; Leininger, T.; Malrieu, J.-P. J. Chem. Phys. 2001, 114, 10252;(b) Angeli, C.; Cimiraglia, R.; Malrieu, J. P. Chem. Phys. Lett. 2001, 350, 297; (c) Angeli, C.; Cimiraglia, R.; Malrieu, J. P. J. Chem. Phys. 2002, 117, 9138. (31) Balint-Kurti, G. G.; Martson, C. C. J. Chem. Phys. 1989, 91, 3571. (32) Hatcher, E.; Soudackov, A.; Hammes-Schiffer, S. J. Phys. Chem. B 2005, 109, 18565. (33) In the case of both 1 and 2, as the HO2− molecule diffuses away from the Ni(III) center, the O22− π*α spin drops in energy relative to the Ni(3dπ) orbital becoming quasidegenerate with the O22− π*β spin, which places the SOMO exclusively on the Ni(III) center. (34) At the B2PLYP level of theory, the following values are obtained: VET = 7920 cm−1, Tna = 6.33 cm−1, Tad = 31.7 cm−1, and Tsc = 21.2 cm−1. (35) Under the conditions used for the stopped-flow experiments examining O2− disproportionation kinetics, a kinetic simulation based on rate constants derived from Scheme 3 indicated that the fourcoordinate species is an important reactant for only the final 2% of the reaction time. (36) Skone, J. H. Quantum Mechanical Methods for Calculating Proton Tunneling Splittings and Proton-Coupled Electron Transfer Vibronic Couplings. Ph.D. Dissertation, The Pennsylvania State University, College Station, PA, 2008. (37) Hazra, A.; Skone, J. H.; Hammes-Schiffer, S. J. Chem. Phys. 2009, 130, 054108.

results for moderately complex systems. This observation is significant because it shows that this treatment can be effectively utilized for medium and large systems not amenable to more costly levels of theory (i.e., double-hybrid DFT and MP2 methods vs CAS-SCF and coupled-cluster methods). Lastly, as previously shown, neither the O2− reduction reaction facilitated by 1 nor 2 could be considered strictly adiabatic or nonadiabatic. Instead, examination of the values obtained for Tsc demonstrated that these reactions are either adiabatic with significant nonadiabatic character or nonadiabatic with significant adiabatic character; however, an examination of the KIE values obtained using Tsc versus Tna and Tad suggests that the tunneling matrix elements at the adiabatic (for 1) or nonadiabatic limits (for 2) could be employed to obtain reasonable estimates of Tsc. Also, when the tunneling matrix element at the “incorrect” limiting case was used, larger deviations from Tsc were observed relative to the “correct” limit. This suggests that for complex systems where the construction of accurate PESs required for the calculation of Tsc would not be feasible, use of tunneling matrix elements at the two limits and subsequent comparison with experimental data can still provide useful insight into the fundamental physics of the system under investigation.

■

ASSOCIATED CONTENT

■

AUTHOR INFORMATION

* Supporting Information S

Cartesian coordinates for the TS structures of 1···O2− and 2··· O2− and plots of the temperature dependence of the KIE values. This material is available free of charge via the Internet at http://pubs.acs.org.

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We thank Prof. Sergey Varganov (UNR) for many insightful discussions and a thorough reading of this manuscript and Prof. Kent Ervin (UNR) for the code to his kinetic analysis software. We also thank the anonymous reviewers for their helpful suggestions. This work was supported by the National Science Foundation (CHE-1362662).

■

REFERENCES

(1) Shearer, J. Acc. Chem. Res. 2014, 47, 902. (2) Broering, E. P.; Truong, P. T.; Gale, E. M. Biochemistry 2013, 52, 4. (3) DuPont, C. L.; Neupane, K.; Shearer, J.; Palenik, B. Environ. Microbiol. 2008, 10, 1831. (4) Bryngelson, P. A., Maroney, M. J. Nickel Superoxide Dismutase. In Metal Ions in Life Sciences; Sigel, A., Sigel, H., Sigel, R. K. O., Eds.; John Wiley & Sons, Ltd.: Chichester, England, 2007; Vol. 2, pp 417. (5) Youn, H.-D.; Kim, E.-J.; Roe, J.-H.; Hah, Y.-C.; Kang, S.-O. Biochem. J. 1996, 318, 889−896. (6) Youn, H.-D.; Youn, H.; Lee, J.-W.; Hah, Y.-C.; Kang, S.-O. Arch. Biochem. Biophys. 1996, 334, 341−348. (7) Lee, J. W.; Roe, J. H.; Kang, S. O. Methods Enzymol. 2002, 349, 90−101. (8) Wuerges, J.; Lee, J.-W.; Yim, Y.-I.; Yim, H.-S.; Kang, S.-O.; Carugo, K. D. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 8569. (9) Barondeau, D. P.; Kassmann, C. J.; Bruns, C. K.; Tainer, J. A.; Getzoff, E. D. Biochemistry 2004, 43, 8038. 5460

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461

Article

The Journal of Physical Chemistry B (38) Hammer, T.; Coutinho-Neto, M. D.; Viel, A.; Manthe, U. J. Chem. Phys. 2009, 131, 224109. (39) Schröder, M.; Gatti, F.; Meyer, H.-D. J. Chem. Phys. 2011, 134, 234307.

5461

DOI: 10.1021/acs.jpcb.5b02640 J. Phys. Chem. B 2015, 119, 5453−5461