Molecular basis of the mechanism of thiol oxidation by hydrogen peroxide in aqueous solution: challenging the SN2 paradigm



通过活性氧物质(reactive oxygen species)氧化含有细胞硫醇的化合物如谷胱甘肽和蛋白质Cys残基,被认为在多种生物过程中起重要作用,包括信号转导,酶活性调节,蛋白质通道,转录因子和抗氧化反应。 过氧化氢(H2O2)在许多细胞类型中产生,作为对各种细胞外刺激的反应,并且可以作为普遍存在的细胞内信使.特别是H2O2对活性蛋白Cys的双电子氧化是关键氧化还原信号和调节过程中的事件.






RS + H2O2 RSO + H2O (eq. 1)

据报道,该反应通过SN2机制进行,但尽管其重要性,但在原子水平上并未完全理解该反应。在这项工作中,我们利用最先进的混合量子经典(QM-MM)分子动力学模拟,阐明了硫醇氧化硫醇对甲硫醇盐模型系统的反应机理。我们的结果表明溶剂在定位反应物中起关键作用,在反应的第一阶段存在显着的电荷再分布,并且在达到过渡态之后在H 2 O 2氧原子之间存在氢转移过程。这些观察结果与该反应的SN2机制假设相矛盾。具体地,我们的结果表明,该反应是由微负电荷的过氧化氧通过对负硫原子的亲电子攻击而在产品中变得甚至更负的趋势驱动的。这与SN2机制不一致,后者预测质子化的亚磺酸和羟基阴离子是稳定的中间体。没有发现这些中间体。相反,反应直接进行到未质子化的次磺酸和水。

Molecular basis of the mechanism of thiol oxidation by hydrogen peroxide in aqueous solution: challenging the SN2 paradigm



The oxidation of cellular thiol containing compounds, such as glutathione and protein Cys residues, is considered to play an important role in many biological processes. Among possible oxidants, hydrogen peroxide (H2O2) is known to be produced in many cell types as a response to a variety of extracellular stimuli and could work as an intracellular messenger. This reaction has been reported to proceed through a SN2 mechanism, but despite its importance, the reaction is not completely understood at the atomic level. In this work we elucidate the reaction mechanism of thiol oxidation by H2O2 for a model methanethiolate system using state of the art hybrid quantum-classical (QM-MM) molecular dynamics simulations. Our results show that the solvent plays a key role in positioning the reactants, that there is a significant charge redistribution in the first stages of the reaction, and that there is a hydrogen transfer process between H2O2 oxygen atoms that occurs after reaching the transition state. These observations contradict the SN2 mechanism hypothesis for this reaction. Specifically, our results indicate that the reaction is driven by a tendency of the slightly charged peroxidatic oxygen to become even more negative in the product via an electrophilic attack on the negative sulfur atom. This is inconsistent with the SN2 mechanism, which predicts a protonated sulfenic acid and hydroxyl anion as stable intermediates. These intermediates are not found. Instead, the reaction proceeds directly to unprotonated sulfenic acid and water.


Keywords: thiol containing compounds, hydrogen peroxide, oxidation, SN2, redox signaling, redox regulation

Go to:


Oxidation of cellular thiol containing compounds such as glutathione and protein Cys residues by reactive oxygen species, is considered to play important roles in a wide array of biological processes including signal transduction, regulation of the activity of enzymes, protein channels, transcription factors and antioxidant responses.1-5 Hydrogen peroxide (H2O2) is produced in many cell types as a response to a variety of extracellular stimuli and could work as a ubiquitous intracellular messenger.6 Particularly, the two-electron oxidation of reactive protein Cys by H2O2 is a key event during redox signaling and regulation.7-9 Though this reaction has been productively studied, the exact mechanism remains poorly understood.


There is consensus about thiolates being much more prone to oxidation than protonated thiols10 and that they react with the protonated form of hydroperoxides.11 The reduction of H2O2 by thiolates in aqueous solution has been described as a bimolecular nucleophilic substitution (SN2) in which the thiolate acts as the nucleophilic moiety generating the breaking of the O–O peroxide bond and displacing the HO– leaving group.11 However, recent theoretical studies have demonstrated that the product of the reaction would be H2O instead of OH–,12-15, which is inconsistent with a SN2 mechanism. In this context, Chu and Trout have investigated the reaction between dimethylthiol and H2O2 via DFT calculations, analyzing the effect of the aqueous environment including up to three explicit water molecules.13 In successive works, using an IMOMO methodology, Cardey and Enescu modeled the oxidation of methanethiolate (CH3S–) and Cys-thiolate (Cys-S–) by H2O2, 14,15 including solvent effects through a Polarizable Continuum Model (PCM).


The overall oxidation reaction has been reported to occur through:

RS + H2O2 RSO + H2O (eq. 1)

We explore the reaction by means of QM-MM MD simulations using an umbrella sampling scheme.26 This approach allowed us to obtain thermodynamical information such as the free energy profile in



Where KaRSH is the acidity constant of the thiol and KaH2O2 is the acidity constant of H2O2, which has a pKa value of ~11.7,16 and therefore is fully protonated at near physiological pHs. addition to microscopic insight about electronic structure changes throughout the reaction. The results presented in this work suggest that the process cannot be considered as a simple substitution since the product formation involves a hydrogen atom transfer.


The data reported herein provides a detailed microscopic view of the thiolate oxidation reaction by H2O2 in water, and is the first step for future studies related to the mechanisms of thiolate oxidation in different protein environments.


Go to:


Initial Survey of the System in Vacuum and Implicit Solvent


In order to obtain valid starting structures for further analysis, to obtain information about the energy surface and the mechanism of the reaction under investigation, and to carry out a methodology evaluation, we performed several electronic structure calculations using Gaussian 03.27 The structures of reactants complex (RC) (H2O2/CH3S–), products complex (PC) (H2O/CH3SO–) and transition state (TS) were optimized both in vacuo and in aqueous solution via the PCM approach,28 at different levels of theory: HF, PBE, B3LYP, MP2, employing a double-zeta plus polarization (dzvp) Gaussian basis set.29 Frequency calculations were performed in all cases. Aiming to investigate if one or more water molecules could be involved in the reaction mechanism, we also performed IRC calculations at the PBE/dzvp level of theory including one and four water molecules in the QM system, both in vacuo and with the PCM method.


The geometry of RC obtained with the B3LYP functional was used as the starting point for a first approximation of the energy profile using the PCM solvation model, in which the reaction coordinate was taken to be the distance between the sulfur atom (S) and the reactive oxygen in the peroxide (Or) in order to obtain partial charges necessary to perform the classical MD simulations which are required to equilibrate the systems, as described below.


QM-MM molecular dynamics simulations


The actual QM-MM simulations were carried out using our own developed program (for details on the QM-MM scheme see refs.30,31). The solute was embedded in a 24 Å cubic box, containing 490 explicit water molecules simulated using the TIP4P mean-field potential.32 The constraints associated with the intermolecular distances in the solvent were treated using the SHAKE algorithm.33 The Lennard-Jones parameters (ε and σ) for the quantum subsystem atoms were 0.250, 0.172, 0.170 and 0.016 kcal/mol, and 3.563, 3.399, 3.250 and 1.069 Å, for S, O, C and H, respectively. For the QM region (CH3S– and H2O2 atoms), computations were performed at the generalized gradient approximation (GGA) level, using the PBE combination of exchange and correlation functionals, with a dzvp basis set for the expansion of the one-electron orbitals. The electronic density was also expanded in an auxiliary basis set and the coefficients for the fitting were computed by minimizing the error in the Coulomb repulsion energy.29


All the QM-MM MD simulations were run for at least 4 ps and employed the Verlet algorithm to integrate Newton's equations with a time step of 0.2 fs. Initial configurations were generated from preliminary 100 ps classical equilibration runs in which the solute was treated classically as a rigid moiety with the RESP charges described above. At t=0, the classical solute was replaced by a solute described at the DFT level, according to the methodology described above. An additional 2 ps of equilibration was performed using the QM-MM scheme. During the simulations, the temperature was held constant at 300 K using the Berendsen thermostat.33 The solute and the rest of the system were coupled separately to the temperature bath.


We employed an umbrella sampling scheme, choosing as reaction coordinate the difference between the Or–Ow and the S–Or distances, which was sampled from -1.5 to 1.5 Å. The harmonic potentials used had spring constants of about 3000 N/m. 31 simulation windows were employed to obtain the free energy profiles.


All dynamics visualizations and molecular drawings, were performed with VMD


Go to:


The free energy profile


After Gaussian geometry optimizations and first initial survey of the system, we performed molecular dynamics simulations using the QM-MM scheme as described above. The obtained free energy profile is depicted in Figure 1.


An external file that holds a picture, illustration, etc.

Object name is nihms-358075-f0001.jpg

Figure 1

Free energy profile obtained by umbrella sampling. Free energy (kcal/mol) is plotted versus reaction coordinate (Å). An illustrative picture with the reaction coordinate distance components and atom names as referred in the text is also shown.


The free energy profile obtained in aqueous solution shows a free energy barrier of ~18 kcal/mol, which is in good agreement with previous experimental18 and theoretical13-15 data available for low molecular weight thiol oxidations by hydrogen peroxide. As expected, the reaction is clearly exergonic, with a change in free reaction energy of about -52 kcal/mol, also in agreement with previously reported information.14


In order to get an estimation of the reliability of DFT in predicting the energetics, we performed characterizations of RC, TS, PC, using different methods (HF, PBE, B3LYP and MP2), both in vacuo and in aqueous solution via the PCM approach, with a dzvp basis set (Table 1 summarizes the results). Frequencies calculations were carried out in all cases.




Reaction Energetics (kcal/mol) for isolated and PCM solvated species


Method  Eact ΔGact     Ereact     ΔGreact

HF/dzvp (in vacuo)       50.2 50.1 -37.1      -36.1

HF/dzvp (PCM)     43.8 41.9 -43.8      -44.2

PBE/dzvp (in vacuo)     15.1 13.9 -33.1      -32.5

PBE/dzvp (PCM)    10.5 8.4   -39.0      -39.7

B3LYP/dzvp (in vacuo) 19.8 18.0 -35.3      -35.3

B3LYP/dzvp (PCM)       15.1 13.3 -42.3      -38.7

MP2/dzvp (in vacuo)    22.7 21.0 -43.1      -43.1

MP2/dzvp (PCM)   14.7 14.8 -52.0      -51.8

It is worth noticing that all structures optimized in this section are very comparable within the different methodologies applied (see Supplementary Information Table 1). As expected, the use of PBE leads to underestimation of barriers in comparison with those obtained with B3LYP or MP2.35


We found that the inclusion of the environment by means of PCM has a significant effect within all models tested here, lowering about 20% the computed values. Furthermore, for all the calculations including the environment effect, except for HF calculations (overestimated barriers), activation barriers are in the range of experimental reports for this kind of reactions.18


However, this apparent qualitative agreement in the computed energetics may be misleading, because both the isolated species and the PCM solvated species exhibit geometrical features in the RC structure that are completely unrealistic for describing the actual situation in aqueous solution. Specifically, the inspection of the geometries reached with both strategies show that, while the vacuum and PCM structures of TS and PC are qualitatively similar to those sampled in the QM-MM simulations, the RC adopts significantly different conformations when explicit waters are present (see Supplementary Information Table 1). With little differences, every geometry optimization performed for the isolated and PCM solvated RC, resulted in a “symmetric” conformation, where the distances between the S atom and both O atoms of the peroxide are practically the same; the same fact is observed for the distances between the S atoms and the H atoms of the peroxide, resulting in a cis-H2O2 (Figure 2a). In contrast, the average conformation of the RC sampled by our QM-MM scheme is an “asymmetric” one, where the Or atom is much closer to the S atom than Ow, and the peroxide adopts a trans configuration (Figure 2b). Consistent with the work described by Chu et al.,13 in the presence of 2-3 explicit water molecules, specific solute-water interactions can form favoring the “asymmetric” arrangement. Since neither in vacuo nor implicit solvents environments allowed these specific interactions, the “symmetric” conformation in which solute atom contacts are maximized is the result of such calculations. Consequently, the RC energies achieved with these two different approaches are not comparable.


An external file that holds a picture, illustration, etc.

Object name is nihms-358075-f0002.jpg

Figure 2

RC conformations comparison. a) PBE/dzvp/PCM geometry of the RC. Note that both peroxide O atoms and their respective H atoms are in a quasi equidistant position from the S atom, in a “symmetric” conformation. b) average structure of the RC sampled with our QM-MM scheme. Oxygen atoms from the peroxide are not equidistant from the S atom, the angle between S-Or-Ow is ~50° larger than the angle obtained through in vacuo and PCM optimizations and the H2O2 shows a trans configuration.


In order to explore the possibility of water molecules being involved in the reaction mechanism, we performed IRC calculations including 1, and 4 water molecules, in the QM system, both in vacuo and with the PCM model. The presence of waters molecules modifies the geometry of the RC in a similar way of that observed in the QM-MM simulation. However, as reported by Chu et al. in a similar model system,13 no significant differences in the reaction mechanism and energetics are observed by including water molecules (see Table S2 and Figure S1 of Supporting Information).


Reaction evolution


Three major events occur during the reaction: the breaking of O–O bond of the peroxide, the formation of the bond between the sulfur (S) and the reactive oxygen (Or) and the transfer of the hydrogen (Ht) bound to Or to the other oxygen of H2O2 (Ow). Our QM-MM scheme allows us to get a microscopic insight into the reaction evolution and mechanism (see the 3D animation available in Supporting Information for illustration).


A close examination of the early TS structure reveals an alignment between S, Or and Ow (see Supplementary Information Table 1). Additionally, Ht is still attached to Or, as shown in the distance versus reaction coordinate plot (Figure 3). This analysis also indicates that the transfer of Ht takes place “downhill”, barrierless, after the system has reached the TS. Therefore, even though the product sulfenic acid appears in its unprotonated form, from a kinetic point of view the leaving group of the reaction could be considered as a hydroxyl, which is in agreement with reaction rates correlation with the conjugated acid of the leaving group pKa.36


An external file that holds a picture, illustration, etc.

Object name is nihms-358075-f0003.jpg

Figure 3

Bond length evolution during the reaction. The distances S–Or , Or–Ow, Or–Ht, Ow–Ht (Å) as a function of the reaction coordinate are depicted using black, red, blue, and light green lines, respectively. The TS region (as determined in Figure 1) is indicated by a dark green box.


The dependence of the Mulliken population over different moieties upon reaction is shown in Figure 4. The system negative charge is in the RC, localized mostly in the S atom, while oxygen atoms of H2O2 possess negative charges of roughly -0.5 e each. Based on experimental results which indicate that all low molecular weight thiols exhibit similar pH independent oxidation rate constants, it has been suggested that charge transfer is not significantly involved in the formation of the transition state.19 However, our results indicate a significant charge rearrangement during the first steps of the reaction resulting in negative charge from the RC distributed across both peroxide O atoms and the S atom. This does not contradict experimental rate constants, since decreases in the solvation in the TS (see below) can cause a significant decrease in the observed dependence of the rate on the basicity of the nucleophile.37 Although in the PC, both Or and Ow atoms possess a charge of about -1 e as would be the case in a homolytic rupture, the reorganization of oxygen charges is not entirely symmetric. Ow suffers a faster decrease (reaching a final status consistent with a water oxygen atom charge) while Or does not reach its final charge until it is bonded with S, acquiring most of the negative charge of the CH3SO– moiety. The methyl group preserves its close to zero charge along the entire process.


An external file that holds a picture, illustration, etc.

Object name is nihms-358075-f0004.jpg

Figure 4

Charge evolution during the reaction. Mulliken charges (e) of S atom (yellow line), Or atom (red line), Ow atom (blue line) and methyl group (black line) are plotted versus the reaction coordinate (Å). The TS region (as determined in Figure 1) is indicated by a dark green box.




Previous studies have modeled solvation effects in this kind of reaction using either continuum models or the inclusion of a small number of explicit water molecules by means of geometry optimizations.13-15,38 Our QM-MM MD approach allows us to obtain a more realistic picture of solvation in bulk water at room temperature. We have analyzed the solvation structure along the process, to monitor how solvation patterns may affect the reaction energetics. Radial correlation functions of selected atoms with water oxygen atoms from the RC (left panel), TS (middle panel) and PC (right panel), along with representative snapshots of each step, are shown in Figure 5.


An external file that holds a picture, illustration, etc.

Object name is nihms-358075-f0005.jpg

Figure 5

Solvation structure. Up: radial correlations functions of S (black line), Or (red line) and Ow (blue line) with water oxygen atoms from the RC (left panel), TS (middle panel) and PC (right panel). Down: representative snapshots of the solvation structure of the RC (left panel), TS (middle panel) and PC (right panel).


Radial correlation functions of S atom show that this moiety losses its hydrophilicity along the reaction. The opposite effect is observed in both oxygen atoms. While in the RC the oxygens are poorly solvated, they become better solvated in the TS due to the early charges redistribution and strongly solvated in the PC. Moreover, the peaks of both peroxide oxygen atoms observed in the TS and PC are located at ~2.8 Å (radial distance), while the S peak in the RC is ~3.2 Å, confirming that it is much more favorable to solvate O– like atoms than the S– one. These observable facts lead to a PC much better solvated than the RC and the TS.




We present here an integrated QM-MM approach for the oxidation of CH3S– by H2O2, which allow us to get microscopic dynamical information of this reaction in aqueous environment. The energetics obtained are in good agreement with previous experimental18 and theoretical13-15 reported data and also with quantum calculations presented herein.


We found that the transfer of the hydrogen atom from the reactive peroxidatic oxygen to the non-reactive one takes place after reaching the TS and without an energy barrier; thus, a water molecule appears in the products, suggesting that in aqueous solution the protonated form of the product sulfenic acid is never produced. Although the calculations performed in this work suggests that no water molecules are chemically involved in the process, the solvent plays an important role in the reaction, particularly in the reactants structure. Furthermore, the activation energy is due to two main events: the alignment of the sulfur atom with the peroxidatic oxygen atoms and a significant charge redistribution, which correlates with important changes in the solvation profile. This phenomenon emphasizes the key role of the environment in modulating thiol reactivity, which can be explicitly evaluated by our QM-MM scheme.


In summary, our microscopic insight is not in agreement with the SN2 paradigm because the product formation involves a hydrogen atom transfer which precludes the possibility of a simple substitution. Even though the initial attack of the reactant thiolate to an oxygen atom of H2O2 is an attack between two negative charged centers, the peroxidatic oxygen electrophilicity is the true driving force in this process as evidenced by significant electron donation from the thiolate to the peroxide in the TS.


Thiol groups in cysteine containing proteins usually react with peroxides with similar rate constants as low-molecular weight thiols.39 However, there are some highly reactive thiols which reduce peroxides surprisingly fast. This is the case for the reaction of hydrogen peroxide reduction by thiol-dependent peroxidases such as peroxiredoxins. Fast reactivity was initially considered to be related to the low peroxidatic thiol pKa that assured thiolate availability at physiological pHs. However, it is now widely recognized that increased thiolate availability cannot explain differences in reactivities by factors of ~106-107 as have been reported for H2O2 reduction by peroxidatic thiols in some peroxiredoxins compared with the uncatalyzed reaction. Very recently, active site microenviromental factors leading to transition state stabilization have been considered.9,19,23 The results presented here emphasize the importance for the protein environment that surrounds these thiols, to contribute in the alignment of the S atom with the peroxide O atoms in the TS,23,40,41 so that a substantial charge reorganization can take place. Moreover, differences involving aqueous and specialized protein environments could produce major changes in the reaction mechanism, perturbing the energetics of PC, TS and RC, and/or modifying the reactants electrophilicity. For example, in the case of peroxiredoxins it has been proposed that the presence of nearby residues such as arginine, proline and threonine, may form an active hydrogen network that could assist in placing the substrate and making a simple substitution mechanism possible.23 Information about the atomistic detailed mechanism of thiol containing compounds oxidation with other oxidants such as peroxynitrite is still scarce. In this context, experimental and theoretical studies are underway in our laboratories, to shed light on the mechanisms of thiol oxidation and the broad modulation of thiol reactivity in proteins towards oxidants and electrophiles.



Chem Res Toxicol. Author manuscript; available in PMC 2013 Mar 19.

Ari Zeida,a Ryan Babbush,a,b Mariano C. González Lebrero,c Madia Trujillo,d Rafael Radi,d and Darío A. Estrina

a Departamento de Química Inorgánica, Analítica y Química-Física and INQUIMAE-CONICET, Facultad de Ciencias Exactas y Naturales, Universidad de Buenos Aires, Ciudad Universitaria, Pab. 2, C1428EHA Buenos Aires, Argentina

c IQUIFIB-Dpto. Química Biológica, Facultad de Farmacia y Bioquímica, Universidad de Buenos Aires, Buenos Aires, Argentina

d Departamento de Bioquímica and Center for Free Radical and Biomedical Research, Facultad de Medicina, Universidad de la República, Av. Gral Flores 2125, CP 11800, Montevideo, Uruguay.

bPresent address: Department of Chemistry, Harvard University, 12 Oxford Street, Cambrigde, MA, USA.


Molecular basis of the mechanism of thiol oxidation by hydrogen peroxide in aqueous solution: challenging the SN2 paradigm  https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3308630/