- Arrhenius Laboratory, Department of Organic Chemistry, Stockholm University, Stockholm, Sweden
When DFT has been applied to study mechanisms of redox processes a common procedure has been to study the results for many different functionals. For redox reactions involving the first row transition metals, this approach has given very different results for different functionals. The conclusion has been that DFT cannot be used for these reactions. In the meantime, results with strong predictability have been generated, most noteworthy for photosystem II, where all DFT predictions have been verified by experiments performed later. In order to obtain these predictive results using DFT, an alternative, systematic approach has been used, where the key differences between the results for different functionals can be rationalized by using a single parameter, rather than using the very large number of differences in the functionals.
I. Introduction
The results using density functional theory for molecules, have been continuously criticized the past decades. In particular, the results using b3lyp (Becke, 1993) have been shown to have severe errors in several cases. For main group elements, most of these cases were found to be corrected by the addition of dispersion corrections (Schwabe and Grimme, 2007). The cases where dispersion was needed could relatively easily be identified, and the majority of calculations made using b3lyp were not affected by dispersion in a significant way. However, molecules containing transition metals have always been considered particularly difficult to treat, and the use of DFT has continued to be criticized in these cases. DFT using different functionals were shown to give very different results for many redox reactions, in particular for molecules containing first row transition metals. The selection of different functionals for different reactions has become popular, but this approach obviously suffers from user bias, and lack of predictability.
When different functionals are compared, it is important to do this in a systematic way. During the past decade it has been possible to identify the most sensitive parameter for redox enzymes, and this is the amount of exact exchange. By identifying this single parameter, it is possible to get away from the non-systematic approach of testing a large number of functionals, which differ in many irregular ways. The best procedure found is to start with the usual b3lyp functional, which has 20% exact exchange, and decrease this percentage in steps.
A breakthrough for that approach came in the studies of photosystem II (Siegbahn, 2009, 2013). The catalyst is an Mn4Ca complex, which shows most of the typical problems of using DFT. For example, the Mn(III) to Mn(IV) redox energy is very sensitive to the fraction of exact exchange used in hybrid methods like b3lyp. For every change of one percent of the exact exchange, the redox energy changes by one kcal/mol (Siegbahn and Blomberg, 2014). This means that for b3lyp using 10% exact exchange, the redox potential may change by 10 kcal/mol compared to the one using the usual 20%. Since there are three oxidations of Mn(III) in the catalytic cycle the difference becomes 30 kcal/mol from the start to the end. Still, it was possible to make new and certain predictions for water oxidation in PSII, that were confirmed in detail by experiments performed later.
II. Computational Details
The general computational approach used here has been applied and evaluated for a large number of metallo-enzymes (Blomberg et al., 2014). Density functional calculations are performed on large cluster models of the active sites (170–300 atoms). The active sites in the enzymes discussed here have two or more transition metal ions, and the models include these atoms, plus the first and sometimes also the second shell ligands. As mentioned in the introduction, the b3lyp functional (Becke, 1993) with a varying fraction of exact exchange is used in the present study. A strong argument for using b3lyp as a reference functional, is that this is the hybrid DFT functional with the smallest number of parameters. In fact, the sensitivity of the results for b3lyp is essentially dependent on only one parameter, the amount of exact exchange. Geometries are optimized using a double zeta basis with polarization functions on all second row atoms, and with a few atoms close to the truncations fixed to the X-ray coordinates. The reason not to use a larger basis set for the geometry optimization is based on a large amount of experience gathered during the past three decades, see for example Siegbahn (2001). In fact, an even smaller basis set would probably be accurate enough. More accurate energies for the optimized structures are obtained from single point calculations using a larger basis set, the lacv3p+ basis for the metal ions (Jaguar, 2009), and the large cc-pvtz(-f) basis set for the rest of the atoms. In the recommended use of b3lyp, the fraction of exact exchange has generally been set to 15% (Reiher et al., 2001) in previous studies. Empirical dispersion corrections according to Grimme (Schwabe and Grimme, 2007; Grimme et al., 2010), and solvent effects from the surrounding protein using the self-consistent reaction field (SCRF) approach are included in the energetic results described below. When H2 is released and N2 becomes bound in nitrogenase, the gain (loss) of translational entropy of about 10 kcal/mol is also included. The inclusion of entropy is very important in these cases but not elsewhere. Zero-point corrections are taken from the Hessians, calculated at the same level as the geometry optimizations. For nitrogenase, these effects were taken from smaller models. The Jaguar 7.9 program (Jaguar, 2009) is used for the nitrogenase geometry optimization and for all the calculations with the larger basis set, and the Gaussian 09 package (Gaussian, 2010) is used for the cytochrome c oxidase geometry optimizations and for the Hessian calculations.
The computational procedure described above has been kept as similar as possible to what has been used the past decades. The reason is that a large experience has been gained of this approach, and this knowledge is very useful when the accuracy of the predictions obtained is judged.
III. Nitrogenase
Nitrogenase is the main enzyme in nature that catalyzes nitrogen reduction from the air. The core of its catalytic cofactor is shown Figure 1 (Kim and Rees, 1992; Spatzal et al., 2011). It has seven irons and one molybdenum connected by sulfide bridges. Species containing vanadium and with all-iron exist. A redox potential of –1.6 V, the lowest in nature, is used for the reduction. Quite surprisingly, the oxidation state of the cofactor is not particularly low, which would normally have been expected for a strongly reducing complex. Instead there are four Fe(III) and three Fe(II). Molybdenum is in the Mo(III) state.
Almost all experimental information concerns the ground state before reduction. However, there is one notable exception which is a key for finding the mechanism by DFT model calculations. It has been shown by EPR that the reduced state that activates N2, termed E4, has two bridging hydrides (Hoffman et al., 2013). These two hydrides were found to be removed as H2 in a reductive elimination process, which is directly followed by the binding of N2. The process was found to be easily reversible by changing the pressure of hydrogen and nitrogen. This means that the states involved must be nearly isoenergetic.
The experimental suggestion for the structure of E4 the central carbon remains unprotonated and has two bridging hydrides (Lukoyanov et al., 2016). Since E4 was suggested to appear after four reductions of the experimentally characterized ground state in the catalytic cycle, see Figure 1, there should be two remaining protonations. They were suggested to be on the sulfides. This meant that the redox state for E4 should be the same as for the ground state with four Fe(III), a very surprising situation.
The results obtained for nitrogenase after four reductions of the ground state in Figure 1 are shown in Table 1. The corresponding structures are shown in Figures 2, 3. For each functional, the results for six different structures are listed. These structures were taken from the best ones obtained previously, where several hundred structures were compared. They were obtained with 20% exact exchange at the lacvp* level, and were then compared with a large cc-pvtz(-f) basis set at the 15% level. For each functional the geometries were optimized, which turned out to be important. The first structure, termed C,H−, has an unprotonated carbon and one hydride. There are three protonations of the sulfides. In the second structure, termed C,2H−, one of the protons on the sulfides has moved to a hydride position, so there are two hydrides and two protonated sulfides. This is the one of the six structures that corresponds to the experimentally suggested structure. In the third structure, termed CH, the hydride in the first structure has moved to carbon, and there are still three protonated sulfides remaining. The fourth structure, termed CH2,H−, has a doubly protonated carbon, one hydride and one protonated sulfide. In the fifth structure, the hydride has moved to carbon to form a CH3 group, and there is one protonated sulfide. Finally, in the sixth state, the two hydrides in the second structure has been removed to form a free H2 molecule.
Table 1. Relative energies (kcal/mol) for states obtained after four reductions from the ground state of nitrogenase, using density functionals with different fractions of exact exchange.
The main conclusion that can be drawn from the results in Table 1 is that the structure, termed C,2H−, does not fit the experimental observations. This is true for all functionals in the table. Already at this point, the important conclusion can be made, that the E4 state does not look as suggested experimentally. However, the by far strongest argument against the experimental suggestion comes from the results for the sixth structure. Experimentally, it is known that two hydrides should be removed as H2 essentially thermoneutrally. As can be seen in the table, the removal of H2 is much too exergonic for all functionals to fit that observation. With 20%, the removal is exergonic by 54.6 (22.6 + 32.0) kcal/mol for 15% exergonic by 48.0 (14.7 + 33.3) kcal/mol, for 10% exergonic by 43.3 kcal/mol and for 0 % exergonic by 34.0 kcal/mol.
It is also found that all functionals prefer a protonated carbon, suggested in previous DFT studies (Siegbahn, 2016). For 20, 15, and 10%, this preference is very large. For 0%, the preference is less pronounced, but it is still there by 11.4 (23.4–12.0) kcal/mol.
It should be added that results with three hydrides are not included in Table 1, but have been done. The reason three hydrides are not in the table is that they can be excluded as possible E4 states already from the beginning. It is chemically unreasonable that the first four reductions of the cofactor should lead to an oxidation of the metals rather than a reduction, in particular, since the lowest redox potential of –1.6 V is used. The already high oxidation state of the ground state, with four Fe(III), would be six Fe(III) with three hydrides. The possibility that there should be five Fe(III) and a Mo(IV) is equally unreasonable chemically. Furthermore, for 20, 15 and 10% exchange the structure with three hydrides is very high in energy compared to the other structures. Even for the functional with 0%, the best solution obtained is higher in energy than the CH2H− structure by +6.1 kcal/mol. However, that value is somewhat uncertain since the spin-coupling with three hydrides is very different from the other cases, and includes several low-spin-coupled Fe(III). This leads to a very large number of possibilities for the spin-state. Five spin-couplings were tried, all higher in energy than the CH2H− structure. They were selected to have the sums of the spins on Fe2-Fe4 opposite to the ones for Fe5-Fe7. Two of them were the ones optimal for the ground state and for the reduced states found previously (Siegbahn, 2016). The last argument against structures with three terminal hydrides is that they are incompatible with the EPR experiments (Hoffman et al., 2013), which clearly show two bridging hydrides.
There are other points of interest in Table 1. For 15 and 20% exchange, the structure with two hydrides are quite high in energy compared to the first structure, by +22.6 and by +14.7 kcal/mol, respectively. Furthermore, they are also terminal, not bridging, hydrides in contradiction to experiments. The energy of this structure is somewhat lower for 10%, but still +7.2 kcal/mol higher than for the starting structure. For 0%, the situation is different. Here the structure with two hydrides is lower than for the first structure by -12.0 kcal/mol, but the two hydrides can still be removed to form H2 with a gain of 34.0 kcal/mol, as mentioned above. The lowest state for 20 and 15% is the CH3 structure, but for 10 and 0% it is the CH2,H− structure, with a margin of –1.0 kcal/mol for 10 % and with –5.3 kcal/mol for 0%.
In summary, there are a few clear conclusions. First, none of the functionals prefers an unprotonated carbon. Also, all functionals bind the hydrides very poorly. To remove the hydrides as H2 is strongly exergonic in all cases, by more than 30 kcal/mol, when experiments show that this should be almost thermoneutral. This means that the experimentally suggested structure for E4, does not agree with the results for any functional and can be ruled out by a large margin by DFT.
In a very recent paper, a new mechanism for H2 release and N2 binding in nitrogenase has been suggested (Raugei et al., 2018). A non-hybrid method was used. A mechanism was suggested in which the two hydrides in E4 endergonically form a locally bound H2 molecule. To avoid the problem with the very large computed exergonicity when H2 is released, the key to their mechanism is that this bound H2 molecule could only be released with a very high barrier. If the barrier is lower than 18 kcal/mol, there would be no protonation of N2, but a high barrier should prevent H2 from being released. In their mechanism, the endergonic binding of N2 was found to reduce the barrier for releasing H2 by a significant amount, and H2 can then be released. There are many bound H2 complexes in the literature, but none of them behaves like suggested in Ref. (Raugei et al., 2018). In all the published cases there is at most a weakly bound H2, which can be released without additional barriers apart from the endergonicity. A search for a bound H2 was initiated in the present study using 0% exchange (non-hybrid), and a picture is obtained which is very similar to the ones previously published in the literature. A weak local minimum for a bound H2 is obtained and releasing it from that minimum goes over at most a very small barrier, less than 5 kcal/mol. The release is quite exergonic including a gain of entropy of about 10 kcal/mol.
Since DFT rules out an unprotonated carbon structure, another structure has to be found for E4. If no other structure can be found, the conclusion must be that no version of DFT can handle nitrogenase, a very unlikely scenario. It has been suggested in previous DFT studies that the lowest energy structure in Table 1, is merely a starting structure for catalysis (Siegbahn, 2016). The first four reductions of the ground state structure in Figure 1 would then be just an initial activation process, done only once, before catalysis starts leading to a new E0 state. Following that suggestion, four additional reductions would lead to the actual E4 structures. The results for the E4 structures, determined in this way, are shown in Table 2 for the same four functionals as discussed above. The results are also displayed in Figure 4. The results for six structures are shown, following the conclusions of the EPR experiments. The first entry, termed 2H−, has two hydrides, a CH3 ligand and three protonated sulfurs, altogether eight protonations of the cofactor. The second entry, termed “H-H re TS” shows the barrier for the re mechanism (reductive elimination) of two hydrides to form H2. The third entry, termed “H-H hp TS,” shows the barrier for the hp (hydride, proton) mechanism, where one hydride and one proton form H2. In the fourth entry, termed “- ,” the two hydrides have been removed to form a free hydrogen molecule. The fifth entry, termed “- ”, differs from the second one by a rotation of the homocitrate ligand. This rotation was found in earlier studies (Siegbahn, 2018) to be required for binding N2. In the final sixth entry, termed “+ N2–,” N2 binds to the fifth structure.
Table 2. Relative energies (kcal/mol) for states obtained after eight reductions from the ground state of nitrogenase, using density functionals with different fractions of exact exchange.
Figure 4. The results for the theoretically suggested E4 structures are shown for different fractions of exact exchange used , see also Table 2.
The results for the functional with 20% agree very well with the EPR experiments. H2 can be removed with a very small energy gain of –2.1 kcal/mol. The barrier for H2 elimination using the re mechanism is 13.1 kcal/mol, well within the range required by experiments. Even more importantly, the barrier using the hp mechanism of 17.3 kcal/mol is 4.2 kcal/mol higher than for the re mechanism, leading to the required N2 activation rather than H2 production. The following homocitrate rotation is only uphill by +1.4 kcal/mol. In the next step, N2 binds exergonically to the third structure by –0.1 kcal/mol. This means that the release of H2 and the binding of N2 should be easily reversible as observed experimentally. There is a minor discrepancy, since there should be some driving force (exergonicity) for this process, but the error is small for such a complicated process. Going to 15%, the discrepancy is somewhat larger. The release of H2 is exergonic by -6.6 kcal/mol and N2 binding is endergonic by +4.3 kcal/mol. However, the discrepancy to experiments is not alarming. Again, the barrier eliminating H2 by the re mechanism is preferred over the one using the hp mechanism, now by 6.2 kcal/mol. For the functional with 10%, the discrepancy to experiments is significantly increased. Most noteworthy, the binding of N2 is now endergonic by + 11.2 (9.3 + 1.9) kcal/mol. However, the barrier for the re mechanism is still lower than the one for the hp mechanism. Finally, for the functional with 0%, the discrepancy to experiments increases further. For example, the binding of N2 is now endergonic by +14.2 kcal/mol. Furthermore, the barrier for the hp mechanism is now lower than for the re mechanism, leading to production of H2 rather than protonation of N2.
IV. Cytochrome C Oxidase
The membrane bound enzyme cytochrome c oxidase (CcO) catalyzes the reduction of molecular oxygen to water as the last step in the respiratory chain in aerobic organisms. The chemistry occurs in an active site consisting of a high-spin heme group, a copper complex, labeled CuB, and a redox active tyrosine, referred to as the binuclear center (BNC). The electrons are delivered to the BNC from a reduced cytochrome c, located on the positive side of the membrane. The protons are transferred to the BNC from bulk water on the negative side of the membrane. Molecular oxygen binds to the reduced BNC, with heme-Fe(II), CuB(I) and TyrOH. The O-O bond is cleaved in a single reaction step, yielding a four electron oxidized BNC, with heme-Fe(IV)=O, CuB(II)OH and a neutral TyrO-radical (Proshlyakov et al., 1998; Fabian et al., 1999). The rest of the catalytic cycle consists of four reduction steps, each taking up one electron and one proton to the BNC, leading back to the reduced state with two new water molecules. The overall energetics of the reduction process is obtained from the difference in reduction potential between the electron donor and the acceptor, molecular oxygen. With cytochrome c as the electron donor, reduction of one oxygen molecule becomes exergonic by 50.7 kcal/mol (2.2 eV) (Brzezinski, 2004). A significant part of this free energy is conserved as an electrochemical gradient across the membrane, which in turn is used by another enzyme, ATP-synthase, to produce ATP, the energy currency of the cells. Two processes contribute to the gradient buildup, one is the electrogenic chemistry (taking the electrons and the protons from opposite sides of the membrane), and the other is the so called proton pumping, which means that the chemistry is coupled to proton transfer across the entire membrane. The largest group of CcO:s are known to pump one proton per electron, i.e., four protons per oxygen molecule (Brzezinski and Gennis, 2008; Kaila et al., 2010). The mechanism for the proton pumping, i.e., how to couple the transfer of one electron to the active site with the uptake of two protons from the negative side of the membrane, is still under debate (Rich, 2017).
The process of oxygen reduction in CcO has been studied in detail using density functional theory. A mechanism for the O-O bond cleavage step was suggested based on computational results at an early stage (Blomberg and Siegbahn, 2006, 2010), and it was later confirmed in a combined experimental and computational study (Poiana et al., 2017). Another result from the computational studies concerns the mechanism of proton pumping, for which it has been suggested that the redox active tyrosine in the active site plays an essential role (Blomberg, 2016). Experimental support for the suggested pumping mechanism is that an active site tyrosine is conserved in all families of CcO's (Hemp et al., 2006). To understand the proton pumping in CcO it is essential to know the energetics of the individual reduction steps in the catalytic cycle, which depends on the reduction potential of the active site cofactor that is reduced in each particular step. Experimental investigations have indicated that the four different active site reduction potentials involved are quite different, and that only two of them seem to be large enough (about 0.8 V) to afford proton pumping (Wikström and Morgan, 1992; Kaila et al., 2010). Also the BNC potentials have been studied computationally, with result that partly differ from the experimental measurements (Blomberg and Siegbahn, 2015a), see further below. In all these computational studies the b3lyp type of functional was used as described above in Computational details.
To further test the reliability of the results obtained for CcO, it was decided to systematically investigate the calculated proton coupled reduction potentials for the active site cofactors by varying the fraction of exact exchange in the b3lyp potential, using the model shown in Figure 5. It is noted that each of the proton coupled reduction potentials corresponds to the formation of a new O-H bond in the active site. The individual reduction potentials can therefore be estimated from the strengths of the different O-H bonds. The strength of an O-H bond is not sensitive to a distant surrounding, since the charge is not changed by the addition of a (H+,e−)-couple, which means that reasonably sized models (150–200 atoms) can be used in the calculations. For each functional (different percentage of exact exchange) new structures are optimized, and since the most recently used procedure for CcO has been to optimize the structures with dispersion, the D3 dispersion correction with parameters from the original b3lyp-D3 functional (Grimme et al., 2010) was used in both geometry optimizations and energy calculations. Since results for the functional with 15% were also obtained using structures optimized with 20% (not reported) it could be noted that the final energetics differs only slightly (<2 kcal/mol) from the ones reported below.
Figure 5. The CcO model used in the calculations, showing the active site cofactors in the four electron oxidized state: CuB(II)-OH, TyrO(288)-radical and heme a3-Fe(IV)=O.
To estimate the energetics of the reduction steps of the CcO catalytic cycle, the sum of the energy of an electron transferred from cytochrome c plus a proton transferred from bulk water is needed. This energy is parameterized to reproduce the experimental results for the overall reaction (50.7 kcal/mol, see above). By subtracting the parameterized energy from each O-H bond strength, the exergonicity relative to the cytochrome c donor is obtained, and by comparing to the midpoint potential of cytochrome c (0.25 V) the midpoint potential can be estimated for each cofactor. The first result to be discussed concerns the proton coupled reduction potential for CuB, for which the experimental measurements have given low values, 0.2–0.4 V (Jancura et al., 2006; Brand et al., 2007; Vilhjámsdóttir et al., 2018). In contrast, the previous calculations indicate a much higher value (0.9–1.0 V) during catalytic turnover (Blomberg and Siegbahn, 2015a), and an explanation for the low experimental values has been suggested based on the computational results (Blomberg and Siegbahn, 2015b). As can be seen in Table 3, all functionals support a large reduction potential for CuB, with results that vary only slightly and in a somewhat irregular way with the fraction of exact exchange in the functional, 0.97-1.19 V. Interestingly, using a small model complex it was also shown that the b3lyp* (15%) results for the O-H bond strength in CuB(I)OH2 agrees to within one kcal/mol with CCSD(T) results (Blomberg and Siegbahn, 2015a).
Table 3. Calculated energetics as a function of the amount of exact exchange for the reduction of CuB (taking up one electron and one proton) in the catalytic cycle of cytochrome c oxidase.
As mentioned above, experiments have shown that two of the active site reduction potentials are about 0.8 V (Kaila et al., 2010), which is large enough for proton pumping in the corresponding reduction steps. These two steps concern reduction of the tyrosyl radical and heme-Fe(IV)=O, respectively (Kaila et al., 2010). The calculated results for these two cofactors using the different functionals are reported in Table 4. For the tyrosyl radical the experimental value is 0.82 V, and all functionals with a fraction exact exchange (10–20%) are in reasonable agreement with experiment, 0.73–0.86 V. Again, the reduction potential does not vary very much with the exact exchange, and even the functional without exact exchange (0%) is not too far off, with a value of 0.60V. For the heme-Fe(IV)=O reduction the situation is quite different. The experimental value is 0.76 V (Kaila et al., 2010), and here the functional with 15% gives the best agreement, with a value of 0.68 V. The result with 20%, 1.01 V, is also in reasonable agreement with experiment. The result without exact exchange (0%) disagrees qualitatively with experiment, the calculated value is –0.22 V, and also with only 10% exchange the calculated value, 0.45 V, is quite far from the experimental value. An important factor for the strong variation of this reduction potential is that there is a change of spin coupling in the heme-Fe(IV)=O → heme-Fe(III)OH transition, where iron is low-spin coupled in heme-Fe(IV)=O and high-spin coupled in heme-Fe(III)OH.
Table 4. Calculated energetics as a function of the amount of exact exchange for two of the four reduction steps (taking up one electron and one proton) in the catalytic cycle of cytochrome c oxidase.
V. Conclusions
By using the present systematic approach of DFT, a few major conclusions can be drawn for the nitrogenase mechanism. In Table 1, the results are given for different E4 structures obtained after four reductions of the ground state, as has been suggested experimentally. The most important result here is that no functional gives results consistent with the conclusions of the experimental EPR study (Hoffman et al., 2013). In particular all functionals give a very large exergonicity for releasing H2. This means that the experimentally suggested structure can not be supported by any DFT functional. The discrepancy to experiments is very large, independent of fraction of exact exchange used. Even the non-hybrid method (0%), gives very poorly bound hydrides, which can be removed with a gain of 34 kcal/mol, when experiments indicate that the release should be nearly thermoneutral. In contrast, the theoretically suggested structure for E4, obtained after eight reductions of the ground state, shows a much more reasonable agreement with the experimental EPR information, in particular, with 20 and 15% exchange. Since DFT indeed found a structure that agrees well with what is known for the E4 state, a conclusion that DFT totally fails for nitrogenase is very farfetched. If that would have been the case, it would be the first example among the many redox reactions studied, that would show that behavior (Blomberg et al., 2014). If no structure would have been found that agreed with experimental information for E4, the conclusion could have been different. A very recent study with many different functionals did not find any functional that gave a preference for the experimental structure (Cao et al., 2018), in agreement with the prediction made here.
For CcO, the calculations show that the active site cofactor CuB has a large midpoint potential, about 1 V, a result obtained regardless of the fraction of exact exchange in the functional. This is in contrast to the much lower potentials obtained in experimental measurements, but in agreement with experimental observations on proton pumping. For another active site cofactor, heme-Fe(IV)=O, the functional with 15% exact exchange gives the best agreement with experimental observations for the midpoint potential, and the functional without exact exchange gives qualitatively wrong results.
In summary, the results for the different functionals show that the best agreement with experiments, for both nitrogenase and cytochrome c oxidase, is obtained with 15–20% exact exchange in the functional. 10% exact exchange gives results slightly worse, and the use of a non-hybrid functional (0% exact exchange) gives qualitatively wrong potential surfaces in both cases. These results are in line with experience gathered during the past two decades (Blomberg et al., 2014). This experience can be used in two ways. It could either be used to calibrate the exact exchange fraction based on some well-known experimental fact for that system. Or, the difference between the results for 15 and 20% could be used as an estimate of the error in the calculations. The latter approach has mainly been used in our previous studies.
Author Contributions
PS wrote the paper and did nitrogenase calculations. MB wrote the paper and did cytochrome c oxidase calculations.
Funding
This work was supported by the Swedish Research Council (grant numbers 2015-04104 and 2016-03721). Computer time was provided by the Swedish National Infrastructure for Computing.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
References
Becke, A. D. (1993). Density-functional thermochemistry. III. The role of exact exchange. J. Chem. Phys. 98, 5648–5652. doi: 10.1063/1.464913
Blomberg, M. R. A., and Siegbahn, P. E. M. (2015a) How cytochrome c oxidase can pump four protons per oxygen molecule at high electrochemical gradient Biochim. Biophys. Acta. 1847, 364–376. doi: 10.1016/j.bbabio.2014.12.005
Blomberg, M. R. A. (2016). The mechanism of oxygen reduction in cytochrome c oxidase and the role of the active site tyrosine. Biochemistry 55, 489–500. doi: 10.1021/acs.biochem.5b01205
Blomberg, M. R. A., Borowski, T., Himo, F., Liao, R.-Z., and Siegbahn, P. E. M. (2014). Quantum chemical studies of mechanisms for metalloenzymes. Chem. Rev. 114, 3601–3658. doi: 10.1021/cr400388t
Blomberg, M. R. A., and Siegbahn, P. E. M. (2006). Quantum chemistry applied to the mechanisms of transition metal containing enzymes-cytochrome c oxidase a particularly challenging case. J. Comp. Chem. 27, 1373–1384. doi: 10.1002/jcc.20448
Blomberg, M. R. A., and Siegbahn, P. E. M. (2010). Quantum chemistry as a tool in bioenergetics. Biochim. Biophys. Acta 1797, 129–142. doi: 10.1016/j.bbabio.2009.10.004
Blomberg, M. R. A., and Siegbahn, P. E. M. (2015b). Protonation of the binuclear active site in cytochrome c oxidase decreases the reduction potential of CuB. Biochim. Biophys. Acta. 1847, 1173–1180. doi: 10.1016/j.bbabio.2015.06.008
Brand, S. E., Rajagukguk, S., Ganesan, K., Geren, L., Fabian, M., Han, D., et al. (2007). F. A new ruthenium complex to study single-electron reduction of the pulsed OH state of detergent-solubilized cytochrome oxidase. Biochemistry 46, 14610–14618. doi: 10.1021/bi701424d
Brzezinski, P. (2004). Redox-driven membrane-bound proton pumps. Trends Biochem. Sci. 29, 380–387. doi: 10.1016/j.tibs.2004.05.008
Brzezinski, P., and Gennis, R. B. (2008). Cytochrome c oxidase: exciting progress and remaining mysteries. J. Bioenerg. Biomembr. 40, 521–531. doi: 10.1007/s10863-008-9181-7
Cao, L., Caldararu, O., and Ryde, U. (2018). Protonation and reduction of the FeMo cluster in nitrogenase studied by quantum mechanics/molecular mechanics (QM/MM) calculations. J. Chem. Theory Comput. 14, 6653–6678. doi: 10.1021/acs.jctc.8b00778
Fabian, M., Wong, W. W., Gennis, R. B., and Palmer, G. (1999). Mass spectrometric determination of dioxygen bond splitting in the “peroxy” intermediate of cytochrome c oxidase. Proc. Natl. Acad. Sci. U.S.A. 96, 13114–13117.
Grimme, S., Anthony, J., Ehrlich, S., and Krieg, H. (2010). A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 132, 154104. doi: 10.1063/1.3382344
Hemp, J., Robinson, D. E., Ganesan, K. B., Martinez, T. J., Kelleher, N. L., and Gennis, R. B. (2006). Evolutionary migration of a post-translationally modified active-site residue in the proton-pumping heme-copper oxygen reductases. Biochemistry 45, 15405–15410. doi: 10.1021/bi062026u
Hoffman, B. M., Lukoyanov, D., Dean, D. R., and Seefeldt, L. C. (2013). Nitrogenase: a draft mechanism. Acc. Chem. Res. 46, 587–595. doi: 10.1021/ar300267m
Jancura, D., Berka, V., Antalik, M., Bagelova, J., Gennis, R. B., Palmer, G., et al. (2006). Spectral and kinetic equivalence of oxidized cytochrome c oxidase as isolated and “activated” by reoxidation. J. Biol. Chem. 281, 30319–30325. doi: 10.1074/jbc.M605955200
Kaila, V. R. I., Verkhovsky, M. I., and Wikström, M. (2010). Proton-coupled electron transfer in cytochrome oxidase. Chem. Rev. 110, 7062–7081. doi: 10.1021/cr1002003
Kim, J., and Rees, D. C. (1992). Structural models for the metal centers in the nitrogenase molybdenum-iron protein. Science 257, 1677–1682. doi: 10.1126/science.1529354
Lukoyanov, D., Khadka, N., Dean, D. R., Raugei, S., Seefeldt, L. C., and Hoffman, B. M. (2016). Photoinduced reductive elimination of H2 from the nitrogenase dihydride (Janus) State Involves a FeMo-cofactor-H2 intermediate. Inorg. Chem. 56, 2233–2240. doi: 10.1021/acs.inorgchem.6b02899
Poiana, F., von Ballmoos, C., Gonska, N., Blomberg, M. R. A., Ädelroth, P., and Brzezinski, P. (2017). Splitting of the O-O bond at the heme-copper catalytic site of respiratory oxidases. Sci. Adv. 3:e1700279. doi: 10.1126/sciadv.1700279
Proshlyakov, D. A., Pressler, M. A., and Babcock, G. T. (1998). Dioxygen activation and bond cleavage by mixed-valence cytochrome c oxidase. Proc. Natl. Acad. Sci. U.S.A. 95, 8020–8025. doi: 10.1073/pnas.95.14.8020
Raugei, S., Seefeldt, L. C., and Hoffman, B. M. (2018). Critical computational analysis illuminates the reductive-elimination mechanism that activates nitrogenase for N2 reduction. Proc. Natl. Acad. Sci. U.S.A. 115, E10521–E10530. doi: 10.1073/pnas.1810211115
Reiher, M., Salomon, O., and Hess, B. A. (2001). Reparameterization of hybrid functionals based on energy differences of states of different multiplicity. Theor. Chem. Acc. 107, 48–55. doi: 10.1007/s00214-001-0300-3
Rich, P. R. (2017). Mitochondrial cytochrome c oxidase: catalysis, coupling and controversies. Biochem. Soc. Trans. 45, 813–829. doi: 10.1042/BST20160139
Schwabe, T., and Grimme, S. (2007). Double-hybrid density functionals with long-range dispersion corrections: higher accuracy and extended applicability. Phys. Chem. Chem. Phys. 9, 3397–3406. doi: 10.1039/b704725h
Siegbahn, P. E. M. (2001). Modeling aspects of mechanisms for reactions catalyzed by metalloenzymes. J. Comp. Chem. 22, 1634–1645. doi: 10.1002/jcc.1119
Siegbahn, P. E. M. (2009). Structures and energetics for O2 formation in photosystem II, Acc. Chem. Res., 42, 1871–1880. doi: 10.1021/ar900117k
Siegbahn, P. E. M. (2013). Water oxidation mechanism in photosystem II, including oxidations, proton release pathways, O-O bond formation and O2 release. Biochim. Biophys. Acta 1827, 1003–1019. doi: 10.1016/j.bbabio.2012.10.006
Siegbahn, P. E. M. (2016). Model calculations suggest that the central carbon in the FeMo-cofactor of nitrogenase becomes protonated in the process of nitrogen fixation. J. Am. Chem. Soc. 138, 10485–10495. doi: 10.1021/jacs.6b03846
Siegbahn, P. E. M. (2018). A major structural change of the homocitrate ligand of probable importance for the nitrogenase mechanism, Inorg. Chem. 57, 1090–1095. doi: 10.1021/acs.inorgchem.7b02493
Siegbahn, P. E. M., and Blomberg, M. R. A. (2014). Energy diagrams for water oxidation in photosystem II using different density functionals. J. Chem. Theory Comput. 10, 268–272. doi: 10.1021/ct401039h
Spatzal, T., Aksoyoglu, M., Zhang, L., Andrade, S. L. A., Schleicher, E., Weber, S., et al. (2011). Evidence for interstitial carbon in nitrogenase FeMo cofactor. Science 334, 940. doi: 10.1126/science.1214025
Vilhjámsdóttir, J., Gennis, R. B., and Brzezinski, P. (2018). The electron distribution in the “activated” state of cytochrome c oxidase. Sci. Rep. 8:7502. doi: 10.1038/s41598-018-25779-w
Keywords: density functional theory, redox reactions, nitrogenase, cytochrome c oxidase, exact exchange
Citation: Siegbahn PEM and Blomberg MRA (2018) A Systematic DFT Approach for Studying Mechanisms of Redox Active Enzymes. Front. Chem. 6:644. doi: 10.3389/fchem.2018.00644
Received: 02 November 2018; Accepted: 11 December 2018;
Published: 21 December 2018.
Edited by:
Vicent Moliner, Universitat Jaume I, SpainReviewed by:
Claudio Greco, Università Degli Studi di Milano Bicocca, ItalyKendall N. Houk, UCLA Department of Chemistry & Biochemistry, United States
Copyright © 2018 Siegbahn and Blomberg. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Per E. M. Siegbahn, cGVyLnNpZWdiYWhuQHN1LnNl