Theoretical Study of O2 Reduction and Water Oxidation in Multicopper Oxidases

Recent electrochemical experiments have shown that the reduction of O2 can be driven backward to water oxidation, which is the first case that has been successfully demonstrated for an enzyme. To understand this ability of the enzyme, both the forward reduction and backward oxidation have been studied here. For the forward reaction, a mechanism similar to earlier studies was obtained. All steps of the full catalytic cycle were obtained for the first time, and it was shown that the explicit reduction steps contribute significantly to the rate-limiting step of the O–O bond cleavage. For the backward oxidation reaction, it was found that the mechanism of the O–O bond formation is not just the reverse of the reduction step where the O–O bond is cleaved for a protonated peroxide. The formation of two fully deprotonated oxo groups was found to be important, which leads to a large radical character for one of the oxo groups. For this possibility, it is important that the pKa of the water bound to the cofactor is quite high.


I. INTRODUCTION
Multicopper oxidases are enzymes that reduce O 2 to water using a copper cluster as the active catalyst. 1,2 The model used for the active site is shown in Figure 1. It was built from an X-ray structure obtained at a high resolution of 1.4 Å for CueO (PDB code 1KV7). 3 The active site has three copper centers with mainly histidine ligands, altogether eight of them. There is also a water molecule bound to one of the copper centers, the Cu T2 center. The other two coppers are termed Cu T3 . Electrons are supplied by yet another copper complex, the Cu T1 complex, situated at a distance of about 13 Å from the Cu 3 complex.
The O−O bond cleavage mechanism for multicopper oxidase has been studied extensively by Rulísěk et al., using different theoretical methods such as QM, QM-MM, multireference, and CASPT2 methods, the latter two of ab initio-type. The first of these studies was performed already in 2005 using QM-MM 4 and it suggested that O 2 binds with one of the oxygens bridging between the Cu T2 ion and one of the Cu T3 ions, while the other oxygen binds to the other Cu T3 ion. After the O−O bond cleavage, an O 2− coordinates with all of the three coppers in the center of the cluster. A year later, the same group used multireference ab initio methods to determine the spectroscopic properties of some of the intermediates in the catalytic cycle, obtaining good agreement with the experiments. 5 In 2007, the peroxy adduct was established by a combination of computational and experimental methods. 6 In 2011, another study was performed to determine the transition state for the O−O bond cleavage. The barrier was found to be 60−65 kJ/mol, in good agreement with the experiments. 7 In 2013, a review was published where the findings from theoretical studies of multicopper oxidases were summarized. 8 Finally, in 2015, most of the intermediates in the catalytic cycle were studied by combined quantum-and molecular-mechanical free-energy perturbation methods. 9 A full reaction cycle was proposed. The effect of the interaction between the distant Cu T1 complex and the tri-copper complex was also studied.
An interesting new development for multicopper oxidases appeared last year when it was shown by electrochemistry that the reduction reaction could be reversed to become an oxidation reaction of water, a very important reaction studied intensively at present. This was the first time that the reduction of O 2 could be reversed for an enzyme. A redox potential essentially the same as for Photosystem II of 1.23 V was used, but at a high pH of 10.5. O 2 was demonstrated to be produced, but at a very low rate. An earlier bioelectrochemical study exists but is not as relevant in the present context of studying the reverse reaction. 11 The present study was inspired by the electrochemical experiment. Both the reduction and oxidation reactions are studied using the same methodology as used for many earlier studies of enzyme mechanisms. 12 As part of the efforts to deepen the understanding of enzyme mechanisms, the ability of some enzymes to reverse reactions has been studied. Besides the chemical interest, it has been shown that the calculated ability to reverse reactions is an excellent way to test if a suggested mechanism is plausible or not. It can also be used to test the accuracy obtained by the calculations. The first test done with that purpose was a study of nickel-containing CO dehydrogenase (Ni-CODH), 13 which is known to be driven backward toward the reduction of CO 2 by a minor decrease of the redox potential. 14 Both the water oxidation and reduction of CO 2 are important reactions in the context of the green-house effect. Recently, also the oxidation of H 2 driven backward to the formation of H 2 was studied theoretically by the same type of theoretical methods. 15 The same reversible reaction is at present being studied also for FeFe hydrogenase. The understanding of H 2 formation is also important in the development of using hydrogen as a fuel in connection to water oxidation.

II. METHODS
The same methods and basis sets were used here as the ones used in recent similar studies 14 and are also quite similar to the ones used in the previous decade on many different redox enzymes. 12 The calculations start out by geometry optimizations using the density functional theory (DFT) functional B3LYP. 16 The basis sets (LACVP*) used in these optimizations are of moderate size but have been shown to be adequate if just the energies are of interest. 17 Backbone atoms were fixed from the Xray structure as described before. 18 For the final energies, a large basis set was used with cc-pvtz(f) for the nonmetal atoms and with LAV3P* for the metals. An important modification of the standard B3LYP method has been used. It has been noted that the energies are only sensitive to the fraction of the exact exchange in the functional. Therefore, this fraction has been suggested to be varied from the standard 20 to 15 and 10%. 19 The D3 dispersion correction 20 was used in the geometry optimizations and for the final energies. Solvation effects were obtained using a Poisson−Boltzmann solver, 21 with a dielectric constant of 4.0. Zero-point effects were taken from calculated Hessians. The present calculations have been done using Jaguar 21 and Gaussian. 22 The experimental driving force for the entire reduction cycle using the redox potential of O 2 to water is −31.2 kcal/mol, using a redox potential of 0.8 eV for the reduction of O 2 at pH = 7 and a redox potential of 0.46 V for the Cu T1 reductant. 23 The direct calculation gives a driving force that does not agree with the experiments. The calculated driving force is therefore adjusted to produce the correct driving force by modifying the cost for each reduction step by the same amount for both of them. The amount is the same since the enzyme surrounding is the same for both of them. Using this procedure, the energy for the uptake of a (H + , e − ) couple from the reductant and water is 387.3 kcal/ mol for 15% exact exchange. To obtain the same driving force, values of 385.8 and 388.6 kcal/mol were used for 10 and 20% exact exchange, respectively. Using an experimental redox potential of 0.46 V for the reductant and an energy of 279.8 kcal/ mol for a proton in water leads to a very similar value of 389.1 kcal/mol. The oxidation of water was performed at pH = 10.5 with a redox potential of 1.23 V for the oxidant, 10 which leads to a total driving force of −65.8 kcal/mol. Fitting to that value, a gain of 411.7 kcal/mol for the donation of a (H + , e − ) couple from the Cu 3 cluster is obtained.
The model used for the Cu 3 active site is shown in Figure 1. The starting state of the cluster is the fully reduced cluster with three Cu(I) states. All eight histidines are kept unprotonated. To stabilize the positively charged cluster, there are two negative aspartates with strong hydrogen bonds to two of the histidines. A water molecule, W1 in Figure 1, weakly bound to Cu T2 is also present. An important additional water molecule, W2 in Figure  1, was found to be bound close to the Cu 3 complex. It is strongly hydrogen-bonded to the other water molecule and to Asp112. During the entire reduction process, the calculated binding energy for W2 is significantly larger than 14 kcal/mol, which is the empirical value for the binding of a water molecule in bulk water. 12

III. RESULTS
Multicopper oxidases are used for the efficient reduction of O 2 using electrons provided by enzymes like laccase and ascorbate oxidase. In the present study, this reduction reaction is therefore presented first, even though the study was initiated by an interest in the reverse oxidation reaction of water. 10 The full reduction cycle was studied, where the electrons were taken from Cu T1 and the protons from water at pH = 7, with a total cost of 387.3 kcal/ mol for one (H + , e − ); see above. The energy diagram for the full reduction cycle has not been computationally studied before, and, as shown below, the inclusion of the explicit reduction steps is important for determining the overall rate. The results are first presented using an exact exchange of 15%. They are then compared to the results using 10 and 20%. The mechanism found is shown in Figure 2.
The first question to be answered concerns the protonation state of the Cu 3 complex. All eight histidines were first found to be neutral, in agreement with the previous findings by Ryde et al. 9 The main question is therefore if there is a water molecule or a hydroxide bound to Cu T2 . In the previous study, a low pK a of only 3 was calculated for this water in the fully reduced state (red), indicating that it should be a hydroxide. A hydroxide was, therefore, kept for all states in the reduction mechanism in the previous study. This is surprising since the distance between Cu T2 and the oxygen is as large as 2.36 Å in the X-ray structure, and also surprising since a water molecule would be expected for a Cu(I) state. In contrast, in the present study, a very large pK a of this water, larger than 15, was obtained for the water in the fully reduced state with three Cu(I) states. In fact, this ligand is found to be a water molecule in all of the states involved in the reduction mechanism. As described below, for the oxidation of The Journal of Physical Chemistry A pubs.acs.org/JPCA Article water, a quite high pK a of this water is actually required for a reasonably low barrier of the rate-limiting step. With a water molecule, neutral histidines, and two negative aspartates, this leads to a total charge of +1 for the model used for the fully reduced state. Since the reduction steps involve the uptake of a (H + , e − ) couple, the charge remains the same in all steps of the reaction.
Starting out from the fully reduced state with three Cu(I), the first step of the reduction process is the binding of O 2 . The calculations give a slight endergonicity of +0.4 kcal/mol for this step. The O 2 structure is similar to the one found by Ryde et al., 6 where one oxygen bridges between Cu T2 and Cu T3 and the other oxygen bridges between the two Cu T3 , see Figure 3. A difference here is that the best bonding found here, by a few kcal/mol, is a mirror image of the structure in the previous study. In the second step of the reduction process, O 2 becomes protonated and the copper cluster is reduced. Also, this step was found to be endergonic, now by +4.7 kcal/mol. At that point, the O−O bond is cleaved with a mechanism very similar to the one given by Rulísěk et al. 7 The transition state (TS) structure is shown in Figure 4. The calculated local barrier for this step is 10.9 kcal/ mol, with a large exergonicity of −26.7 kcal/mol. The structure of the product is shown in Figure 5. However, to get the ratelimiting barrier, the entire reduction process has to be investigated. Already at this stage, the total barrier has to include the two previous steps with an endergonicity of +5.1 kcal/mol and a barrier in total of 10.9 + 5.1 = 16.0 kcal/mol. After the cleavage, the three coppers are in the Cu(II) oxidation state. The spins are in the range 0.5−0.6 and they are antiferromagnetically coupled.
To complete the reduction process, three additional reduction steps have to be considered. First, there is a protonation of the oxo group produced in the O−O cleavage step, combined with a reduction of the complex. This step is exergonic by −6.8 kcal/mol, producing two bound hydroxides. In the next step, one of the hydroxides becomes protonated to    The Journal of Physical Chemistry A pubs.acs.org/JPCA Article form a water molecule, which leaves the complex. The Cu complex becomes simultaneously reduced. An empirical value of 14 kcal/mol is used for the binding of a water molecule in bulk water. The exergonicity is −3.8 kcal/mol. Finally, to get back to the fully reduced structure, there is an endergonic step by +1.0 kcal/mol, in which the second water molecule is produced and the copper complex is reduced. Since the final step is endergonic, the total rate-limiting barrier has to include that step also, and it therefore is 16.0 + 1.0 = 17.0 kcal/mol. The energetics for the entire process is shown in Figure 6.
To get an estimate of the accuracy of the calculations, the energies were recomputed with 10 and 20% exact exchange. The results are shown in Figure 7. With 10%, the barrier for the O−O cleavage step goes down from +10.9 to +9.3 kcal/mol. Also, the cost for reaching the reactant for the O−O cleavage goes down from +6.1 to +4.7 kcal/mol, see the figure. This means that the rate-limiting barrier is 14.0 kcal/mol, compared to 17.0 kcal/mol for 15%. As expected, the results for 20% go in the opposite direction. The O−O cleavage barrier is +12.2 kcal/mol and the cost for reaching the reactant is +8.4 kcal/mol. The total ratelimiting barrier thus is +20.6 kcal/mol. The mechanisms are the same for the different percentages. An estimate of the error of the rate-limiting barrier using 15% should be about ±3 kcal/mol, which is fairly normal. The barrier of +20.6 kcal/mol for 20% is somewhat too high compared to the experiments, which indicate a barrier of about 14 kcal/mol. 23 The best agreement is probably found using 10%, but the results for 15% are also quite satisfactory, considering that two reduction steps contribute to the barrier height.
In the previous study by Srnec et al., 7 a barrier of 14−16 kcal/ mol was obtained using B3LYP with 20% exchange. This is somewhat different from the present value of +20.6 kcal/mol. However, in the previous study, the reduction steps and the binding energy of O 2 were not included in the barrier. Using 20% in the present study, these steps amount to +8.4 kcal/mol.
The main objective of the present study was to investigate both the forward reduction reaction and the reverse oxidation reaction of water, which has recently been studied by electrochemistry. 10 To drive the reaction backward, a much higher redox potential of 1.23 V had to obviously be used compared to 0.46 V for the forward reduction reaction. However, this increase was not enough, and the pH also had to be increased from 7.0 to 10.5, which makes it correspondingly easier to release the protons to water. This means that the driving force for the reaction changes dramatically from 31.2 kcal/mol for the forward reduction reaction to 66.4 kcal/mol in the opposite direction, a change of 95.6 kcal/mol. It is remarkable that the enzyme can survive this large change.
In previous studies, for Ni-CODH and the hydrogenases, the mechanisms for the reverse reaction were found to be simply step by step going backward compared to the forward reaction. From the inspection of the reduction diagram in Figure 6, it is clear that this would not work in the present case. The O−O bond-forming step, from the O,OH reactant to the O 2 H product, is uphill by +26.7 kcal/mol. Going over the TS for this step, the barrier is +37.6 kcal/mol, which is prohibitive. It should be emphasized that the barrier for this O−O bond-forming step is dependent neither on the redox potential nor on the pH, and will not be lowered by the changes made. The three oxidation steps, where a (H + , e − ) couple is released, are very exergonic with the high driving force, but this is of no use for overcoming the high barrier.
Another mechanism is therefore required to understand the experimental electrochemical experiment, which undoubtedly produces O 2 . Many unsuccessful attempts were made, but finally a plausible mechanism was found. The first obvious attempts focused on an initial oxidation, which would remove a proton from the bound water and oxidize Cu T2 , but this change did not lower the barrier. Instead, the oxidation reaction was found to start by three oxidation steps, with a water molecule bound to Cu T2 , which are the reverse of the last three reduction steps. It was confirmed that the pK a of the water is high enough to not be  The Journal of Physical Chemistry A pubs.acs.org/JPCA Article deprotonated, even with the much higher redox potential and pH. When the O,OH step is reached, the barrier for the O−O bond formation at this stage is +37.6 kcal/mol, as mentioned above. The mechanism finally obtained is shown in Figure 8.
A fourth oxidation is therefore required, which introduces some Cu(III) character in the Cu-cluster. The lowest energy for this oxidized state is found by removing the proton on the water on Cu T2 (see Figure 9). At this stage, a redox potential of 1.23 V and a pH of 10.5 are possible. In the three initial oxidation steps, it was always found to be much more preferable to remove a proton from one of the added waters. The fourth oxidation step is exergonic by −6.1 kcal/mol. Rather surprisingly, the barrier for the O−O bond formation at this oxidized state, with Cu(III) at Cu T2 , is still much too high. The calculated barrier is 34.4 kcal/mol. The TS structure is very similar to the one for the O− O bond cleavage in the reduction mechanism (see Figure 4).
A solution was found by removing the proton on the bridging OH group, rather than one on the water on Cu T2 . The resulting state has two bridging unprotonated oxo groups and a water molecule on Cu T2 , which otherwise has the same structure as in Figure 9. This O,O state is +6.9 kcal/mol higher than the O,OH one. From this excited state, after four oxidations, the O−O bond can be formed with a local barrier of +16.5 kcal/mol. The TS is shown in Figure 10. It should be noted that the total ratelimiting barrier should include also the excitation energy of +6.9 kcal/mol to reach the reactant. The barrier then becomes +23.4 kcal/mol. This is a quite high barrier, but is in agreement with the observation that water oxidation was found to be very slow. 10 It is possible that the calculated barrier could be a few kcal/mol too high, which is a common finding with the present methods. The full oxidation diagram is shown in Figure 11. The branching point in the fourth oxidation should be noted and was found to be the key for the mechanism.
The calculated spin populations fully explain why there is a much lower barrier for the O−O bond formation of +16.5 kcal/ mol than for the O−OH bond formation of 34.4 kcal/mol. For the O,OH structure, there is a Cu(III) with zero spin on Cu T2 and no spin on the oxygens, while the O,O structure places the oxidation on the oxo group, which has bonds to all the three copper centers. The spin on this oxo group is −0.68, indicating a large radical character, which makes it much more reactive than the oxygens of the O,OH structure, where the oxygen spins are zero. The O,OH state prefers to oxidize Cu T2 since there is an additional negative ligand (OH) on this center.
The high pK a of the water on Cu T2 is the key to the possibility of reversing the mechanism of water oxidation. If the pK a of the water had been much lower, it would be easier to remove the proton to form a negative hydroxide, which would stabilize the Cu(III) character on Cu T2 and would avoid the formation of a radical character on the oxygens, required for the low barrier for the O−O bond formation. Therefore, the pK a enters directly into the stability of the unreactive O,OH compared to the reactive O,O state.    Figure 10 is similar to the one found for biomimetic di-copper complexes. 24 A difference is that in Cu oxidase, a third copper also affects the barrier. It is not surprising that it is more difficult to form the O−O bond in Cu oxidase since the third copper forms a bond with one of the oxygens holding the oxygens apart. It is more surprising that the O−O bond cleavage barrier is also higher for Cu oxidase, 25 even when the mechanism appears quite similar to the one found for the di-copper complexes. It is therefore clear that also the different ligands and the redox properties of the complexes have an effect on the barrier height.
With the TS for the O−O bond formation found (see Figure  10), it is interesting to go back to the reduction cycle to see how this type of TS compares with the one for the O−OH bond cleavage (see Figure 4). Going backward over the TS in the diagram in Figure 11, a barrier of 66.4−46.3 = 20.1 kcal/mol is obtained. This is slightly higher than the O−OH cleavage barrier of 17.0 kcal/mol in Figure 6. Therefore, the conclusion that reduction and oxidation go over different transition states still holds.

IV. DISCUSSIONS AND CONCLUSIONS
The mechanism for the reduction of O 2 to water for multicopper oxidase has been studied, also including the explicit reduction steps. In recent electrochemical experiments, 10 it has been shown that the reduction can actually be reversed to the oxidation of water to form O 2 . This is the first example of an enzyme that has this ability. The reverse oxidation reaction has therefore also been studied here.
For the reduction reaction, a mechanism similar to the one suggested by Ryde and Rulísěk. 8 has been obtained. However, an important difference is that the water ligand on Cu T2 has here been found to have a very high pK a , while in the previous study it was found to be very low, suggesting it to be a hydroxide. Another new finding is that two endergonic reduction steps and the binding of oxygen contribute to the rate-limiting barrier. While the local barrier for cleaving the O−O bond is 10.9 kcal/ mol, the total rate-limiting barrier is 17.0 kcal/mol. The calculated barrier is about 3 kcal/mol higher than the experimental suggestion, which is a usual discrepancy when using the present methodology. An error bar for the calculated barrier has been estimated by varying the amount of exact exchange from 15 to 10 and 20%. It has previously been found that the fraction of exact exchange is the most sensitive parameter in the calculations. 19 The optimal mechanism is the same for these three different fractions.
To experimentally be able to reverse the reaction to oxidation of water, a high redox potential of 1.23 V and a pH of 10.5 had to be used. With these values, a mechanism for the oxidation reaction has here been found. Perhaps surprisingly, it was found that the oxidation mechanism is somewhat different from the reverse one of reduction. While in the reduction mechanism, the O−O bond is cleaved for a protonated peroxide, in the reverse oxidation reaction, the O−O bond is formed between two unprotonated oxo groups. The reason for this difference is that the spins of the two unprotonated oxo groups indicate one of them to have a substantial radical character, which makes it very reactive. In contrast, if the O−O bond is formed between an oxo and a hydroxide ligand, they are both found to have zero spins, with a resulting barrier for the O−O bond formation, which is much higher. To be able to use the formation of O 2 from the two oxo groups, a very high pK a of the water on Cu T2 is necessary. In contrast, in the previous study, a very low pK a was found, which would not have allowed the present oxidation mechanism. It is interesting to note that, given the very different cell potentials and pH, there is in fact no reason why microscopic reversibility should occur.
With the multicopper oxidase reaction, one more case of a redox enzyme mechanism has been shown to be well determined by the present DFT methodology. 19   The Journal of Physical Chemistry A pubs.acs.org/JPCA Article