The influence of protein electrostatics on potential inversion in flavoproteins

Abstract

Biology uses relatively few electron-transfer cofactors, tuning their potentials, electronic couplings, and reorganization energies to carry out the required chemistry. It is remarkable that the potential ordering of two-electron transfer active flavins can be normal (first oxidation at low potential and second oxidation at high potential) or inverted, and the gap between the potentials can be as large as one volt. Analysis based on structural bioinformatics and electrostatics indicates that the ordering of the flavin redox potential is influenced by protein electrostatics. In all 36 flavoproteins examined, the introduction of a negative charge near the flavin in silico increases the extent of potential inversion (by lowering the electrochemical potential of the second electron-transfer step); the introduction of a positive charge near the flavin favors normally ordered potentials. We also find that the addition of positive charges increases the electrochemical potential for the naturally occurring one-electron transition in flavodoxins (between deprotonated hydroquinone and neutral semiquinone) and also increases the second one-electron transition in bifurcating flavins (between anionic semiquinone and fully oxidized flavin). Finally, we find that proximity of a proton acceptor, notably conserved arginine, supports proton-coupled electron transfer because it may act as a proton acceptor, promoting potential inversion. This key arginine residue may enable two-electron transfer chemistry by promoting the proton-coupled electron transfer process over the pure electron transfer process, suggesting how a protein's flavin environment may influence one- or two-electron chemistry in flavoproteins.

---

1 Introduction

Electron transfer (ET) reactions underlie bioenergetics and biocatalysis.^1–5 Flavins are redox cofactors that can perform one- or two-electron chemistry.^6,7 In flavodoxins, the redox cofactor is the flavin mononucleotide (FMN). Flavodoxins carry out one-electron chemistry and have normally ordered potentials (Fig. 1A). That is, the first electron leaving the fully reduced flavin is more strongly reducing than the second electron, although a second redox step is not observed under physiological conditions.⁶ Flavoproteins that carry out two-electron chemistry on spatially separated redox chains are known as electron bifurcating flavins. The potentials of the exiting electrons are often substantially different. In bifurcating flavoproteins, the redox cofactor is flavin adenine dinucleotide (FAD). Fig. 1B illustrates the typical redox states accessible to flavins. In flavodoxins, the FADH⁻ (deprotonated hydroquinone species) species can be oxidized to form FADH (the neutral semiquinone). In bifurcating flavins, FADH⁻ can lose 2 electrons and 1 proton in two steps to form FAD. Bifurcating flavin cofactors are found to have inverted potentials: they release their second electron at a potential that is much more reducing (lower) than the first^6,8 (Fig. 1A). Flavin-based electron bifurcation (FBEB), first reported in 2008,^9,10 has generated intense interest.^{6–8,11–15} The influence of the flavin potentials on electron bifurcation, especially on their possible role in minimizing short circuiting between the low- and high-potential pathways, is of particular interest.^16–18 Inverted one-electron potentials of flavins can differ from the mean (midpoint) potential by as much as 0.5 V. Although potential inversion has been studied in small molecules,^19–23 the molecular origins of potential inversion in flavoproteins are poorly understood. There is no known correlation between the type of flavin performing single or multi-electron transfer and the relative values of the two half-potentials for the E₁ and E₂ steps. In molecules such as dinitrobenzenes, potential inversion was attributed to redox coupled structural rearrangement¹⁹ and hydrogen bonding effects.^16,24 However, these electrochemical studies of potential inversion do not involve flavins.

Fig. 1 (A) Scheme showing the electrochemical potential orderings in normal and inverted regimes. The electrochemical potential of the first electron is higher than that of the second for inverted potentials (i.e., the second electron is more strongly reducing than the first). The electrochemical potential of the first electron is lower than that of the second for normally ordered potentials. E°(1) and E°(2) are the standard electrochemical potentials of the first and second electron steps. (B) Observed redox states for a flavin that can undergo one- and two-electron chemistry. FADH⁻ is the deprotonated hydroquinone species, FADH is the neutral semiquinone species, FAD⁻ is the anionic semiquinone species, and FAD is the fully oxidized species. (C) Isoalloxazine ring of the flavin. Specific R groups define different flavin derivatives (e.g., LFN, FMN, and FAD).

Experimental studies have examined how protein residues near the flavin binding site influence electrochemical potentials.^25–30 Those studies used site-directed mutagenesis to make single residue mutations to alter the hydrogen-bonding interactions near the flavin. The studies found that hydrogen bonding between flavin N5 (Fig. 1C) and the protein strongly influences the stability of the three relevant flavin redox states. O'Farrell et al. found that residue side chain bulkiness destabilizes the flavin semiquinone in all mutants, emphasizing the importance of glycine-61 and flavin N5 for the flavodoxin D. vulgaris.²⁸ For the flavodoxin in A. vinelandii, Alagaratnam et al. found that removing the hydrogen bond between flavin N5 and cysteine (with the C69A mutation) caused the electrochemical potential of the oxidized/semiquinone flavin species to shift to a more positive value and the semiquinone/hydroquinone electrochemical potential to shift to a more negative value.²⁶ These experimental studies do not address the physical origins of potential inversion in flavoproteins. Proton coupling to flavin redox chemistry and the protein environment are more generally thought to influence potential inversion;^31–35 a detailed understanding of how the structure influences flavin potential inversion remains lacking. Earlier studies noted a conserved arginine in bifurcating flavoproteins near N5 of the flavin cofactor (Fig. 1C), and hydrogen bonds to this arginine are thought to provide the source of inverted potentials.^16,36 Earlier studies have also suggested that the flavin redox processes are influenced by the local environment through electrostatic, hydrophobic and π–π stacking interactions.^6,37 Yet, a conclusive link of potential inversion with the protein structure and dynamics is yet to be established. We aim to use electrostatic analysis to explore how the protein electrostatics influence the redox potentials of flavoproteins, including flavodoxins (which perform single-electron chemistry) and bifurcating flavins (which perform two-electron chemistry). By comparing the flavin environments of flavodoxins and bifurcating flavins, we aim to understand if the electrostatic environment could interchange the function of these two kinds of flavoproteins by changing their electrochemical potential ordering. We also investigate the role of the conserved arginine as a possible proton acceptor and the influence of the arginine residue on the flavin potential inversion.

Earlier computational studies explored how protein and cofactor structures influence the electrochemical potentials of redox cofactors. Analysis based on continuum electrostatics,^38,39 hybrid quantum/classical mechanics,^40–51 and ab initio molecular dynamics⁵² was used to study structure–function relations for electrochemical potentials in proteins. Electrochemical potentials of copper proteins,^53–57 iron–sulfur proteins,^58–60 flavoproteins,^61–63 cytochromes,^64–73 heme proteins,^74–77 and proteins with novel cofactors (e.g., the oxygen-evolving complex of photosystem II) were also studied with quantum mechanical/molecular mechanics (QM/MM) and electrostatic methods.^78,79 Theoretical analysis has focused on experimentally measured potentials and on understanding how protein mutations influence the potentials. The electrochemical potentials of azurin and its mutants,^56,57,80,81 rubredoxin,^58,82 lignin peroxidase,⁸³ and T. versicolor laccase T1 (ref. 46) were computed. While computed electrochemical potentials are generally consistent with experimental measurements, fully quantum calculations are beyond reach as they are computationally intensive. Most computational approaches for computing cofactor potentials (and reorganization energies) in proteins use QM/MM methods of MD sampled geometries.⁸⁴ Testing the accuracy of theoretical electrochemical potential calculations in specific proteins requires experimental benchmarks. However, few experimental values of cofactor electrochemical potentials have been reported for electron bifurcating proteins. The protein Nfn1 is the exception.^16,85 The first electron-transfer step in Nfn1 is proton-coupled (PCET), and the second step is simple ET.^16,85 This sequence of reactions is typical of flavin-based electron bifurcation. Earlier studies found that N1 and N5 of the flavin (Fig. 1C) are the sites with the largest change in electron density that accompany flavin oxidation/reduction.^86,87

We denote the electrochemical potentials in the bifurcating proteins as E₁ (for the first PCET electron transfer step between FADH⁻ and FAD⁻) and E₂ (for the second electron step between FAD⁻ and FAD). We refer to the pure one-electron redox step as (between FADH⁻ and FADH), and we define the normal potential regime when E₂ − E₁ > 0. The inverted potential regime is defined when E₂ − E₁ < 0. The electrochemical transitions are observed at the standard electrochemical potentials and and correspond to the following half-reactions, respectively:

The potential difference between the E₂ and E₁ redox steps before mutation is ΔE = E₂ − E₁; after mutation ΔE_mut = E_2mut − E_1mut. The difference between ΔE and ΔE_mut indicates whether the mutation causes the potentials to become more or less inverted: if ΔE > ΔE_mut, the potentials are more inverted upon mutation; if ΔE < ΔE_mut, the potentials are less inverted upon mutation. Potential inversion is calculated for the flavodoxins for the purpose of comparison to bifurcating flavins. The computed two-electron process for the flavodoxins is purely hypothetical. The flavodoxins are one-electron redox proteins.

We find that changing the charges of residues near N1 or N5, or making other chemical changes to the flavin microenvironment near N1 or N5, alters the gap between the first and second electrochemical potentials in both flavodoxins and bifurcating flavoproteins, without inverting the ordering of the two electrochemical potentials. We also find that mutations to the flavin microenvironment can influence the electrochemical potential of the one-electron step ( step) in both flavodoxins and bifurcating flavins. Specifically, our computations find that adding negatively charged residues promotes potential inversion by decreasing the electrochemical potential of the E₂ step (the potential becomes more negative), while the bioinformatics analysis shows that there are more positively charged residues in bifurcating flavoproteins near the flavin binding site. We find that adding positive charges to the flavoprotein environment near the flavin binding site in both flavodoxins and bifurcating flavins increases the electrochemical potential of the E₂ step, thus decreasing the extent of potential inversion, consistent with the simple electrostatics argument on stabilizing or destabilizing reactants and products. We find that the presence of a proton acceptor near the flavin promotes two-electron chemistry and that the presence or absence of a proton acceptor may aid in tuning between one- and two-electron chemistry, provided that the proton acceptor is initially neutral so that it may accept a proton.

2 Computational methods

We used the multiconformation continuum electrostatics (MCCE) method software suite^38,88,89 to calculate flavin electrochemical potentials. We calculated the electrochemical potentials for each flavoprotein in our dataset. Previous studies used MCCE to calculate electrochemical potentials of redox active proteins (see the SI and ref. 38, 88 and 89) by using continuum electrostatics and solving the linear Poisson–Boltzmann equation to describe energetics. The solvent and the protein are treated as continuous dielectric media with charges assigned to the protein atoms from the Amber94 forcefield.⁹⁰ MCCE analysis computes contributions from electrostatic interactions and solvation energy, while also calculating the protonation states of ionizable residues. The protonation states influence the electrostatic interactions within the protein and between the protein and its environment. We used the MCCE approach to calculate electrochemical potentials for FMN in flavodoxins and for the bifurcating FAD in bifurcating flavoproteins. We also used MCCE analysis to study how protein side chain charges influence flavin redox potentials. We compute relative electrochemical potentials for 30 flavodoxins and 6 bifurcating flavoproteins. We used MCCE to compute 3 redox potentials: E₁, E₂, and (see Fig. 1B) for both flavodoxins and bifurcating flavins. E₁ is the PCET step between the deprotonated hydroquinone (DEP) and the anionic semiquinone (ASQ) for both flavodoxins and bifurcating flavins. E₂ is the step linking ASQ to the fully oxidized (OX) species (the second electron transfer step of bifurcating enzymes). is the one-electron step from the DEP species to the neutral semiquinone (NSQ), the one-electron reaction of flavodoxins. Flavin partial atomic charges were assigned using RESP charges⁹¹ calculated using Gaussian 16 (ref. 92) with the B3LYP density functional⁹³ and 6-31+G** basis set.⁹⁴ Atomic charges were generated for each redox species that participates in the electron-transfer half-reaction for each flavoprotein: DEP, NSQ, ASQ, and OX. Reaction ftpl files were created that contain bond information, including bond connectivities and partial atomic charges following the tutorial protocol on the MCCE github page (https://gunnerlab.github.io/Stable-MCCE/) with the pdb2ftpl.py script provided in the /bin directory. Charges were added to each ftpl file from the RESP charges calculated for each flavin cofactor. ftpl files were created for each electron transfer process that we analyzed (E₁, E₂, and ). Lennard-Jones and torsion interaction parameters were taken from the AMBER94 force field;⁹⁰ electrostatic pairwise energies and interactions were calculated by solving the Poisson–Boltzmann equation using the DelPhi software suite.^95–98 For both flavodoxins and bifurcating flavins, we selected an interior protein dielectric constant of 4, as used in many prior electrostatic studies,^{38,73,88,89,99,100} and a dielectric constant of 80 for water. Physical mutations were made in silico by modifying a set of amino acid residues and selecting a side chain rotamer with the highest statistical frequency that avoids steric clashes with the surrounding (frozen) protein (based on analysis with PyMOL^101,102). We performed two different kinds of charge mutations. The first method sampled only the neutral conformer of each residue. This was achieved by modifying the head3.lst file to disable the sampling of charged conformers for residues that are typically polar and charged. The second kind of mutation changes a neutral nonpolar residue to a charged residue. In this case, we change side chain charges without changing the number or placement of atoms. To perform the second kind of charge mutation, the existing MCCE ftpl files for the mutated residues are copied from the MCCE bin directory with existing residue ftpl files and their charges are modified so that the charge (equal to ±1e) is distributed uniformly on all atoms of the residue side chain (see the Results and discussion section). When performing both charge and physical mutations, we assume that the proton lost during the E₁ step exits to solvent.

3 Results and discussion

3.1 Bioinformatic analysis of polarity for residues near flavin binding sites

We first investigated residue compositions near the flavin binding site in flavodoxins and bifurcating flavins. There are thirty flavodoxins and six bifurcating flavoproteins with unique sequences in the protein databank. We analyzed whether or not charged polar residues are found near the flavin (within 12 Å of the isoalloxazine N1 or N5 atoms) in native flavodoxins or FBEB proteins with known structures (Fig. 1C).¹⁰³ A 12 Å radius around N1 or N5 represents the region near the flavin that has been mutated in experimental studies.^26,29 This 12 Å cutoff distance also reflects the weak influence of electrostatic interactions beyond this distance, as the typical interaction energy between two point charges at 12 Å is approximately 301 mV with a dielectric constant of about 4. We first count the number of residues (for all proteins in the dataset) that have heavy atoms within a prescribed distance of N1 or N5. We then divide the count at each distance by the total number of proteins in our dataset (30 for flavodoxins and 6 for bifurcating flavins). This quotient defines an occurrence frequency. The frequency of finding positively or negatively charged residues within 12 Å of N1 or N5 for the flavodoxin and bifurcator sets is shown in Fig. 2. We assumed typical charge states for amino acid side chains:¹⁰⁴ Asp and Glu are polar negative, Arg and Lys are polar positive, and all other residues have zero side chain charge.

Fig. 2 Residue occurrence frequencies for polar residues as a function of N1 and N5 distances in flavodoxins and bifurcating flavins. (A) Positive and negative residues of flavodoxins as a function of distance from N1. (B) Positive vs. negative residue frequency as a function of distance to N5 in flavodoxins. (C) Positive vs. negative residue distances to N1 in the bifurcating flavins. (D) Positive vs. negative residue distances to N5 in bifurcating flavins. (E) Positive vs. negative residue distances to N1 in bifurcating flavins in the absence of arginine being treated as positively charged. (F) Positive vs. negative residue distances to N5 in bifurcating flavins in the absence of arginine being treated as positively charged. We find that the flavin environment in flavodoxins is more negatively charged than the flavin environment in bifurcating flavins, which show more positively charged residues at smaller distances (assuming that arginine is considered positive).

We find that negatively charged residues occur more frequently than positively charged residues at distances from 4 to 12 Å in flavodoxins (see Fig. 2A and B). The closest positive residues to N1 or N5 appear at larger distances of about 9 Å in flavodoxins.

At short distances (R < 6 Å), bifurcating flavins have more positively charged residues near N5 than negatively charged residues (see Fig. 2D). Most bifurcating flavins have a negatively charged residue near N1 (at a distance of 6 Å or less). Negative and positive residues arise with similar frequencies at larger distances (distances greater than about 10 Å from N1 and N5) in bifurcating flavins. At distances of 7–10 Å from N1 and N5, negatively charged residues (e.g., Asp and Glu) occur more frequently than positively charged residues (e.g., Lys and Arg). At distances of 12 Å and more in bifurcating flavins, the frequency of positively and negatively charged residues is the same (whether measuring from N5 or N1; the frequency of finding positive or negative residues at these distances from N5 is within 3% of the N1 values).

The number of negatively charged residues within 7 Å of the N1 site is approximately the same for flavodoxins and bifurcating flavins, suggesting that negatively charged residues may not play a key role in promoting potential inversion. The most significant difference between the cofactor environments in bifurcating flavins compared to flavodoxins is the presence of the positively charged arginine residue near N5 in bifurcating flavins (but not in flavodoxins). This leads us to hypothesize that the presence of positively charged residues could promote inverted potentials, since bifurcating flavins perform two-electron transfer and have inverted potentials.⁶Fig. 2A–D show the positive and negative charged residue counts, assuming standard pK_a values at pH 7. Previous studies suggested that arginine could be a proton acceptor.^16,36 If arginine is protonated under standard conditions (as we assume in Fig. 2A–D), then arginine cannot be a proton acceptor. However, it is well known that amino acid side chain pK_a values in buried positions can vary from their standard canonical values.¹⁰⁵ To account for the possibility that arginine may have a lower pK_a value in bifurcating flavoproteins, we perform the same informatics analysis as in Fig. 2A–D, but now we assume that arginine is neutral (Fig. 2E and F). We find that if arginine is neutral, there are more negative than positive charges found at distances of 5 Å or more, with the largest change in residue occurrence within 7 Å of N5 and 8 Å of N1. From an informatics point of view, there is no significant difference between the bifurcating flavins and flavodoxins as negative charges dominate positive charges from 4 Å to 12 Å, as shown in Fig. 2A, B, E and F. Moreover, there are no proton acceptors in the flavin microenvironment near N1 and N5 in flavodoxins (since they perform one-electron chemistry, the process, proton acceptors are not needed).⁶ As well, there are no negatively charged residues near N5 of the flavodoxins that are within the hydrogen-bonding distance to act as proton acceptors after electron transfer (a proton acceptor is defined as a residue that can accept a proton and adopt a net +1 charge, assuming the canonical amino acid charges noted above). While Alagaratnam et al. studied the hydrogen bond between N5 and a nearby Cys residue, we classify Cys as neutral, because it remains uncharged under standard physiological conditions (pH 7). While Lys and Arg residues must be neutral before accepting a proton, they are positively charged in the informatics analysis indicated in Fig. 2A (green line). We find that up to 8 Å from the flavin, there are no hydrogen bonding residues for the flavin (Lys and Arg) in the flavodoxin environment. As 8 Å is far outside of the hydrogen bonding range, we conclude that flavodoxins have no proton accepting residues near the flavin isoalloxazine ring. The main influence of arginine on the flavoprotein function is to serve as a proton acceptor. We also explored how the first and second redox potentials changed when positive or negative residues were introduced by chemical mutation using MCCE.

3.2 Influence of amino acid residue electrostatics on flavin potential inversion

We studied how the electrostatic environment around the flavin may influence potential inversion. We hypothesize that the proximity of negative charges near the flavin cofactor will enhance the energetic favorability of the one-electron transition (E₂ process). Additionally, we hypothesize that the addition of negative or positive charges is expected to have minimal electrostatic impact on the proton-coupled electron transfer (PCET) step (E₁ process). We first calculated flavin electrochemical potentials for native flavodoxins and bifurcating flavoproteins. Then, we performed two kinds of in silico mutations, as described in the Computational methods section. Charge mutations model the change in the electrostatic environment without changes to the atomistic protein structure, while physical mutations change the atomic structure of the protein (Fig. 3). Details of mutation strategies are described in the SI. We made charge and physical mutations to residues within 8 Å of N1 and N5 in all bifurcating flavins and flavodoxins. An 8 Å cutoff distance from N1 and N5 was chosen for the mutations, as it represents the point at which the most significant differences in residue frequency are observed, as shown in Fig. 2.

Fig. 3 Representation of (A) charge mutations and (B) physical mutations. The purple spheres represent the flavin isoalloxazine ring, the red sphere represents a charged aspartic acid side chain, the green sphere represents a neutral aspartic acid, and the green rectangle represents an alanine residue (neutral charge and different residues compared to the initial aspartic acid residue). In our calculations, we are analyzing how a change in side chain charge near the flavin can tune the electrochemical potential ordering. Changes in side chain charge can be reflected as only changes in charge (charge mutations) or changes in charge and the structure (physical mutations).

The most notable difference between the flavodoxins and bifurcators (from the bioinformatics plots) is that the bifurcating flavins have more positive charges near the flavin compared to the case for the flavodoxins (assuming that arginine is positively charged). However, if arginine is to act as a proton acceptor, it must be neutral. Flavodoxins have more negative than positive charge near N1 and N5 (see Fig. 2A and B). Combining these observations with the fact that bifurcating flavins are known to perform two-electron chemistry and flavodoxins are known to perform one-electron chemistry (and flavodoxins have no arginine residue or other proton acceptors near the flavin cofactor), we hypothesize that the proximity of negative charges near the flavin cofactor will favor the energetics of the one-electron reduction (E₂ process, eqn (2)). This is because the reactant and product differ in net charge by 1 electron in the E₂ electron transfer step, thus increasing the free energy (ΔG⁽⁰⁾) of the reaction. Specifically, adding negative charges to the flavin environment will destabilize the anionic semiquinone reactant due to charge repulsion, contributing to a smaller ΔG⁽⁰⁾ for the reaction. We hypothesize that the addition of negative or positive charges is also expected to have a small electrostatic influence on the proton-coupled electron transfer (PCET) step (E₁ process, eqn (1)), because the reactant and product both maintain a net charge of −1, producing very small changes to the free energy of the reaction. To test the hypothesis that the protein electrostatic environment plays a significant role in modulating the flavin redox potentials, we performed electrostatic analysis for protein mutations, as described below. A summary of the calculations, the question addressed with each calculation, and the results indicating how these mutations influence potential inversion are given in Tables S1–S3 in the SI.

3.2.1 In silico mutations to bifurcating flavoproteins.
3.2.1.1 Charge and physical mutations from nonpolar to polar charged residues near the bifurcating flavin. We first examined how introducing positive and negative charges to the flavin environment may influence potential inversion (charge mutations do not change the atomic structure of the protein). We identified the nonpolar residue nearest to the flavin and selected this residue for charge mutation. We found that introducing negative charges causes the bifurcating flavin potentials to become more inverted (Fig. 4A). Introducing positive charges decreased the potential inversion in all six bifurcating flavoproteins (Fig. 4B). Adding positive charges near the bifurcating flavin reduces the average potential inversion by 10% for charge mutations and 26% for physical mutations (Fig. S5). When adding a negative charge near the bifurcating flavin, the average potential inversion increases by 19% for charge mutations and 34% for physical mutations (Fig. S5). Adding negative charge lowers the electrochemical potential of the E₂ step (the anionic semiquinone to fully oxidized species step) and has minimal impact on the E₁ step, making E₂ − E₁ more negative, thus promoting potential inversion.

Fig. 4 (A) Redox potential orderings after adding negative charges to the bifurcating flavin microenvironment. (B) Redox potential orderings after adding positive charges to the bifurcating flavin microenvironment. The potentials become more inverted as E₂ − E₁ becomes more negative. For all bifurcating flavins, we find that adding negative charges makes the potential ordering more negative (more inverted) while adding positive charges makes the potential ordering more positive (less inverted).

3.2.1.2 Charge and physical mutations from polar charged to nonpolar residues near the bifurcating flavin. We investigated the effect of switching off all negative charges associated with residues within 8 Å of the flavin from the N1 and N5 sites in bifurcating proteins using both charge and physical mutation strategies. We find that neutralizing all protein negative charges causes the E₂ − E₁ value to increase in all bifurcating flavins, as indicated in Fig. 5A. That is, the electrochemical potential ordering becomes less inverted upon removing negative charge, since the flavin electrochemical potentials of E₂ increase upon turning off the protein negative charges. On removing all negative charges within 8 Å of the N1 and N5 sites in the bifurcating flavins, the average percent change for the potential inversion is 22% for the charge mutations and 24% for the physical mutations. The extent of potential inversion is dominated by the change in the E₂ electrochemical potential.

Fig. 5 Redox potential ordering upon mutating negative charges to original PDB structures (unmutated) in the flavin microenvironment in bifurcating flavins (left) and flavodoxins (right). The potentials become more inverted as E₂ − E₁ becomes more negative. For all bifurcating flavins and flavodoxins, we find that removing negative charges from the flavin environment decreases potential inversion.

3.2.2 Mutations to flavodoxins.
3.2.2.1 Charge mutations from nonpolar to charged polar residues in flavodoxins. We next explored how charge mutations that are made to nonpolar residues influence the electrochemical potentials of flavodoxins. Evenly distributing a −1e charge on the side chain across all heavy and light atoms of the nonpolar residue closest to N1 and N5 atoms to the closest residue side chain heavy atom increases the electrochemical potential inversion in all 30 flavodoxins (residues selected for mutation are tabulated in the SI). The same charge mutation (evenly distributing a −1e charge) of a residue near N5 also increases potential inversion in all 30 flavodoxins. Making a positive charge mutation to a nonpolar residue (evenly distributing a +1e charge across all heavy and light atoms on a nonpolar side chain) near N1 or N5 decreases the potential inversion in all 30 flavodoxins, as shown in Fig. 6. Introducing positive charges near the flavin cofactor decreases the average potential inversion by 247%, while introducing negative charge mutations induces an average potential inversion change of 106%.

Fig. 6 (A) Redox potential orderings after adding negative charges near the flavin cofactor in flavodoxins. (B) Potential orderings upon adding positive charges near the flavin cofactor in flavodoxins. The potentials become more inverted as E₂ − E₁ becomes more negative. For all flavodoxins, we find that adding negative charges makes the potential ordering more negative (more inverted) while adding positive charges makes the potential ordering more positive (less inverted).

3.2.2.2 Mutations from charged polar to nonpolar residues near the flavin cofactor in flavodoxins. We mutated all negatively charged residues within 8 Å of N1 and N5 (closest edge-to-edge residue atoms) in all flavodoxins by performing both charge and physical mutations, as with the bifurcating protein structures. We found a similar qualitative result in this case as we did with the bifurcating flavins: neutralizing the negative charges causes the predicted value of E₂ − E₁ to increase in all flavodoxins (i.e., decreasing potential inversion), as shown in Fig. 5B. The E₂ potentials increase when turning off the protein negative charges. By removing all negative charges near the isoalloxazine in flavodoxins, we find that the average percent change for the flavodoxins' potential inversion is −51% for charge mutations; this value decreases by 52% for physical mutations. The change in the E₂ potential dominates the change in extent of potential inversion; shifts in E₁ upon charge mutation are minimal, as found in the bifurcating flavins (see Fig. 4).

3.2.3 Electrostatic influence on the single-electron transfer reaction . We next explored how the flavin electrostatic environment influences the redox step (FADH⁻ → FADH + e⁻). Potential inversion is found in bifurcating flavins, while flavodoxins facilitate one-electron transfer reactions . However, a second electrochemical potential can be measured in flavodoxins, and the two potentials are normally ordered, as illustrated in Fig. 1B. To assess the effects of charge and physical mutations (vide supra) on the first flavin oxidation, we calculated the gap between electrochemical potentials for the step for all flavodoxins and bifurcating flavins. Consistent with our previous analysis of the E₂ step (FAD⁻ → FAD + e⁻), we find that adding negative charges decreases the electrochemical potential for the step, while adding positive charges increases the electrochemical potential, as shown in Fig. 7 and 8.

Fig. 7 Redox potential orderings in flavodoxins upon (A) adding negative charges near the flavin cofactor for the process and upon (B) adding positive charges near the flavin cofactor for the process. For all flavodoxins, we find that adding negative charges makes the electrochemical potential of the process more negative while adding positive charges increases the electrochemical potential of the process.

Fig. 8 Redox potential orderings in bifurcating flavoproteins after (A) adding negative charges near the flavin cofactor for the process and (B) adding positive charges near the flavin cofactor for the process. For all bifurcating flavins, we find that adding negative charges makes the potential ordering more negative (more inverted) while adding positive charges makes the potential ordering more positive (less inverted).

Adding negative charges increases the driving force for one-electron transfer reactions. As discussed above, changes in potential inversion in our mutation studies are attributed mainly to the large differences in the E₂ potentials rather than to shifts to the E₁ potentials. The electrochemical potential changes in both and E₂ follow a similar trend, where adding positive charges decreases the driving force for one-electron chemistry (see the SI for individual charge and physical mutation plots). The electrochemical potential values calculated are relative. In summary, we conclude that negatively charged environments near the flavins promote one-electron chemistry by shifting the potentials to be more negative, favoring the process (FADH⁻ → FADH + e⁻) compared to the corresponding shifts induced in the E₁ potentials process (FADH⁻ → FAD⁻ + H⁺ + e⁻). The bioinformatics analysis indicates that the major difference between the flavodoxin and bifurcating flavin environments is that the bifurcating flavins have more positive charges near the N5 site compared to flavodoxins. However, the positive charges in bifurcating flavins arise from the conserved arginine. This difference in charge suggests that a coarse-grained informatics approach focusing only on the electrostatic environment may be insufficient to determine the origins of potential inversion, as the conserved arginine may play two roles. Arginine may act as a hydrogen bond partner, or may act as a proton acceptor (both with N5). Next, we analyze the role that the conserved arginine may play in tuning flavoproteins to carry out one- and two-electron chemistry.

3.2.4 Summary of charge and physical mutation effects on flavin potential inversion. The computational results described above can be understood in the framework of elementary electrostatics. When the reactant and product have the same charge, a nearby positive charge shifts the energies by about the same amount (note the absence of change in the E₁ values upon changing side chain charges in both heavy and light atoms in SI Fig. S3 and S6). However, when the reactant and product have different charges (as in the E₂ step, where the transition is from an anionic semiquinone state to a fully oxidized species with an accompanying charge change of +1e), the reactant is destabilized by nearby negative charges. For the E₂ step, from anionic semiquinone to fully oxidized species, the product and reactant have different charges, so negative charge near the cofactor destabilizes the reactant (lowering its electrochemical potential). A summary of charge effects on the flavin potentials is summarized in Fig. 9.

Fig. 9 Schematic representation of how adding charges to the flavoprotein environment through both charge and physical mutations alters the flavin electrochemical potentials. (A) Adding positive charge increases the potential of E₂, thus decreasing the potential inversion of the flavin. (B) Adding negative charge decreases the potential of E₂, thus increasing the flavin potential inversion.

Our findings indicate that negative charges near the flavin promote potential inversion, while positive charges near the flavin reduce potential inversion. However, these conclusions are not supported by the bioinformatic analysis shown in Fig. 2. The bioinformatic data indicate the prevalence of positively charged residues in bifurcating flavoproteins near the cofactor (assuming that arginine is positively charged). In contrast, our electrostatic results show that positive charges decrease the driving force for potential inversion. However, the E₁ process in bifurcation is a PCET process, leading us to hypothesize that arginine acts as the proton acceptor, to increase the driving force of the initial step (E₁) of electron bifurcation. By decreasing the positive charge in the system (e.g., by neutralizing a positively charged arginine) near the flavin, we enhance inverted potentials. However, prior studies suggested that the conserved arginine may facilitate two-electron chemistry in bifurcating flavoproteins^16,36 through the flavin forming a hydrogen bond with the conserved arginine. Direct evidence for Arg to act as a proton acceptor has not been reported. We find (Fig. 2) that there are many negatively charged residues within 7 Å of N1 in flavodoxins. This supports the idea that negatively charged residues may play a significant role in modulating one-versus two-electron chemistry. Therefore, we suggest that arginine acts as the proton acceptor and thus promotes two-electron chemistry in electron bifurcation, a hypothesis explored below.

4 Role of arginine in electron bifurcation and potential inversion

4.1 Proton acceptors may promote two-electron chemistry

This section will assess the role of arginine in promoting hydrogen bonding and its potential influence on potential inversion, as well as its ability to facilitate proton dissociation from the flavin, enabling two-electron chemistry by acting as a proton acceptor. We assess the potential roles of arginine assuming that its side chain pK_a is both protonated (standard canonical pK_a) and neutral.¹⁰⁵

4.2 Influence of arginine hydrogen bonding on potential inversion

To evaluate the influence of hydrogen bonding on potential inversion, we compared the relative electrochemical potentials of the wild-type bifurcating flavins with that of arginine in its native orientation to electrochemical potentials obtained after modifying the arginine torsion angle to position its side chain beyond a hydrogen-bonding distance. In the calculations described above, we assumed that the proton leaving the flavin flowed to the solvent. Here, we model arginine as the proton acceptor. We changed the arginine torsion angle to examine the influence of disrupting the hydrogen bonding of the arginine to the flavin on the electrochemical potential of the bifurcating flavin. Our results show that changing the torsion angle of the arginine does not enhance potential inversion, as indicated in Fig. 10A. The main source for altering potential inversion through modifying the arginine dihedral angle is the increase in electrochemical potential of the E₁ step (caused by an increased distance of the proton acceptor from N5), as indicated in Fig. S7. The increased distance between the flavin and arginine would be expected to destabilize the hydrogen bond between the two species. At these larger distances, the protonated Arg is expected to be less stable than it would be at shorter distances, increasing the electrochemical potential of the deprotonated hydroquinone/anionic semiquinone (E₁) process. Modifying the arginine dihedral angle raises the electrochemical potential of E₁ without much influencing the E₂ electrochemical potential; therefore, dihedral angle changes from the wild-type experimentally determined structures promote potential inversion.

Fig. 10 Effect of modifying the conserved arginine dihedral angle on (A) the potential inversion of bifurcating flavin and (B) bifurcating flavin pK_a. We find that modifying the arginine side chain dihedral angle out of the hydrogen bonding range promotes potential inversion, but leads to an increase in pK_a of the bifurcating flavin, both caused by an increase in free energy of the E₁ process reactant (deprotonated hydroquinone).

We hypothesize that arginine dihedral changes from the wild-type structures will increase the pK_a of the bifurcating flavin, hindering proton release from the flavin because of the increased electrochemical potential. MCCE was used previously for pK_a analysis in proteins,^106–111 and we use MCCE to test our arginine dihedral hypothesis by calculating the pK_a values of the bifurcating flavins for the wild-type and modified torsion angles. All pK_a calculations were performed in the same manner as the electrochemical potential calculations described in the Computational methods section, except that in step 4 (Monte Carlo sampling) pK_a values are sampled. These pK_a calculations allow us to assess how changes in hydrogen bonding might influence proton dissociation energy. We find that modifying the arginine torsion angle from the wild-type experimentally determined structures increases the pK_a of the bifurcating flavin (see Fig. 10B). A higher pK_a value of the bifurcating flavin in the modified torsion angle structure (see Fig. 10B) is caused by the loss of stabilizing hydrogen bonds between the flavin N5 and the arginine side chain provided by the conserved arginine in the wild type. This increased pK_a of the modified arginine indicates that the proton is more tightly bound and less likely to dissociate when the arginine side chain is out of the hydrogen bonding range with the flavin N5 atom. This effect translates into a higher free energy requirement for the proton to shift from the flavin to the arginine, suggesting that arginine in the native protein acts as a strong base (since higher K_a corresponds to a lower acid dissociation constant (K_a)), indicating that the proton is less likely to dissociate. Since the free energy change for proton dissociation is given by ΔG⁽⁰⁾ = −RT ln K_a a lower K_a leads to a more positive ΔG⁽⁰⁾, meaning that more energy is required for the proton to transfer. This increased energy barrier confirms that the modified arginine binds the proton more tightly, reinforcing its role as a strong base in the native protein environment as higher pK_a implies lower K_a, indicating a stronger base. These results also suggest a trade-off between hydrogen bonding and potential inversion because increasing the distance between the flavin HN5 proton and the arginine side chain weakens the hydrogen bond, increasing the reactant free energy for the E₁ reaction (deprotonated hydroquinone – neutral arginine), thus decreasing the driving force for proton loss. The weakened hydrogen bonding has a greater influence on the E₁ process than on the E₂ process, decreasing the flavin ET driving force (as indicated by a higher electrochemical potential), while minimally impacting the E₂ electron transfer process (resulting in a lower electrochemical potential for the flavin), as shown in Fig. S7. Arginine atoms selected for torsional modification, original dihedral angles of the wild-type conserved arginine, and modified dihedral angles are detailed in Tables S3 and S4. Changes in the first and second flavin electrochemical potentials upon arginine dihedral angle changes are found in Fig. S7.

4.3 Assessing the likely proton acceptor: arginine versus water

The proton transport pathway from N5 in bifurcating proteins remains unknown. The proton may transfer to nearby water molecules, amino acid residues, or solvent via a residue-mediated pathway, as suggested in earlier studies.¹⁶ To investigate possible proton acceptors coupled to the first electron transfer step and to explore how each possible proton acceptor influences potential inversion, we computed the electrochemical potential of each bifurcating flavin with different corresponding proton acceptors (i.e., water vs. arginine). Since it is well known that amino acid side chain pK_a values can vary widely with protein environments, we investigate the possibility that arginine may serve as a possible proton acceptor, meaning that the arginine will begin in a neutral ionization state before the E₁ step of electron bifurcation can occur. If the proton transfers from the flavin to the arginine, the proton could move to NH1, NH2, or NE (protonated arginine has charge delocalized across its guanidinium group). We described the reactant and product for the E₁ step (deprotonated hydroquinone to anionic semiquinone) with a combined flavin/arginine ftpl file (details of creating this file are found in the SI) to assess which nitrogen of the arginine side chain is the most likely proton acceptor. This analysis allows us to examine how the flavin potential inversion differs with arginine or water serving as the proton acceptor. In this combined file, the reactant has arginine in the neutral form and the product arginine has a proton added to one of the nitrogen atoms of the arginine side chain (we analyzed each of these protonated structures one by one). When we use internal water as the proton acceptor, the conserved arginine is positively charged.

To examine whether or not bulk water can act as the proton acceptor, an explicit water molecule was added to the crystal structure. The oxygen atom of the added water was positioned 3 Å from the proton attached to the N5 atom of the bifurcating flavin, ensuring that there are no steric clashes with other protein atoms, and our calculations find that the computed potentials remained robust with respect to water orientation at this distance. A 3 Å distance was selected as it falls within the typical range for hydrogen bonding. We find that the computed flavin potentials are more inverted when arginine is the proton acceptor compared to water, regardless of which arginine nitrogen accepts the proton (Fig. 11). We find that proton transfer to any arginine side chain nitrogen generates more inverted electrochemical potentials compared to the case of proton transfer to the explicitly added water, as is seen in Fig. 11. The primary evidence for this finding is that the E₁ values for the flavin are lower when arginine is the proton acceptor compared to when an internal water molecule is the acceptor. This enhancement of potential inversion occurs because the larger free energy of the reactant results in a larger driving force for the E₁ process when arginine is the acceptor rather than an internal water molecule. This finding explains the increased potential inversion observed when the proton is transferred to arginine rather than to water, as shown in Fig. 11. Protonating the nitrogen of the conserved arginine closest to N5 produces the greatest inversion compared to the other nitrogens on the arginine side chain or internal water. Although arginine can accept a proton in its neutral state, the kinetics of proton association and dissociation with the flavin are unknown. A protonated arginine would preclude its role as a proton acceptor. If the arginine can release its proton prior to the bifurcation step, it can serve as the proton acceptor. The flavin electrochemical potential calculations had the proton on either the conserved arginine or the added explicit water molecule as the proton acceptors. Fig. 2D shows a preponderance of positive charges near N5 of the bifurcating flavin (based on canonical charges). The pK_a values of protein residues, especially arginine, are known to depend on the protein environment.^105,112,113 Although our calculations suggest that arginine is more likely to serve as the proton acceptor than water (as indicated by the prevalence of flavin inverted potentials), additional experimental data on arginine pK_a values and on flavin proton transfer kinetics would be needed to determine conclusively whether arginine or water is the primary proton acceptor. While we account for both arginine and solvent as potential proton acceptors, the mechanism of proton loss under standard conditions remains unknown.

Fig. 11 Changes in potential inversion for the different possible proton transfer paths in the first electron bifurcation step. The proton transfer paths considered are a water molecule, the closest nitrogen of the conserved arginine near HN5, the 3 nitrogens on the conserved arginine considered one by one and then their potentials averaged as possible acceptors, and water acting as a proton acceptor with arginine removed from the crystal structure. Our calculations indicate that arginine is more likely to act as the proton acceptor rather than water, as potential inversion is more significantly enhanced (potential ordering is more negative) when any nitrogen on the arginine side chain serves as the proton acceptor compared to water.

5 Conclusions

Structural bioinformatics and continuum electrostatics were used to analyze 6 bifurcating flavoproteins and 30 flavodoxins. The analysis indicates how the protein environment influences flavin electrochemical potential ordering. We find that the electrostatic characteristics of the protein are predicted to influence potential inversion significantly, and our findings explain the normally ordered potentials found in flavodoxins compared to the inverted potentials found in bifurcating flavoproteins. Specifically, our results indicate that the introduction of negative charges near the flavin enhances potential inversion by destabilizing the negatively charged anionic semiquinone, while minimally influencing the neutral oxidized flavin, thus promoting potential inversion. Conversely, introducing positive charges near the flavin promotes normally ordered potentials. This relationship between charge distribution and electrochemical potential ordering suggests that the protein environment tunes flavin potentials, favoring either one-electron or two-electron (E₁) chemistry (depending on the presence of a proton acceptor) and influencing potential ordering. Without a proton acceptor, our analysis finds that introducing additional negative charges near the flavin favors one-electron chemistry via the process as the reactant is destabilized (deprotonated hydroquinone). Conversely, in the presence of a proton acceptor, two-electron chemistry is promoted by increasing the driving force of the E₁ process, and adding negative charge near the flavin appears to enhance potential inversion.

We also found that when the conserved arginine of bifurcating flavins acts as a proton acceptor, it contributes significantly to potential inversion (in bifurcating flavoproteins); the potentials are more inverted when arginine, as compared to an explicit water molecule, is the proton acceptor. While we are able to provide design constraints on how to promote flavin potential inversion, we are unable to demonstrate that electrostatics alone accounts for tuning between one-electron chemistry (normally ordered potentials) and two-electron chemistry (inverted potentials). The structural plasticity of proteins, charge mutations, and dynamic interactions between residues may also play roles that are not captured in our analysis. Transient measurement of the protonation state for the conserved arginine in bifurcating flavoproteins and experimental pK_a values would assist in tracking the PCET dynamics.

Although the computations described here provide a framework for understanding electrochemical potential tuning in flavoproteins, comprehensive studies that integrate protein structural flexibility and increased atomistic detail are needed to build further on our findings. However, the current findings offer a starting point for understanding how the physical and chemical properties of proteins drive potential orderings and enable specific catalytic mechanisms in flavoproteins.

Author contributions

NS, PZ, and DNB designed the research; NS performed the research; NS and PZ analyzed the data; and NS, PZ, and DNB wrote the paper.

Conflicts of interest

There are no conflicts to declare.

Data availability

The datasets supporting this article as well as informatics scripts and MCCE inputs can be found on GitHub at https://github.com/nivsingh/electrostatics_code/tree/main/Electrostatics_Data.

Supplementary information: Informatics plots, E₁ and E₂ plots for all mutations, and additional computational details. See DOI: https://doi.org/10.1039/d5sc02960k.

Acknowledgements

This material is based on work supported by the U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences, under award number DE-SC0021144. We graciously acknowledge Dr Junjun Mao and Prof. Marilyn Gunner for very helpful discussions about the MCCE software suite. We thank Dr Kelsey Parker and Dr Jonathon Yuly for helpful discussions.

References

1. H. B. Gray and J. R. Winkler, Proc. Natl. Acad. Sci. U. S. A., 2005, 102, 3534–3539.
2. V. L. Davidson, Acc. Chem. Res., 2008, 41, 730–738.
3. R. A. Marcus and N. Sutin, Biochim. Biophys. Acta, Bioenerg., 1985, 811, 265–322.
4. D. N. Beratan, C. Liu, A. Migliore, N. F. Polizzi, S. S. Skourtis, P. Zhang and Y. Zhang, Acc. Chem. Res., 2015, 48, 474–481.
5. D. N. Beratan, Annu. Rev. Phys. Chem., 2019, 70, 71–97.
6. K. Kayastha, S. Vitt, W. Buckel and U. Ermler, Arch. Biochem. Biophys., 2021, 701, 108796.
7. W. Buckel and R. K. Thauer, Front. Microbiol., 2018, 9, 401.
8. J. L. Yuly, C. E. Lubner, P. Zhang, D. N. Beratan and J. W. Peters, Chem. Commun., 2019, 55, 11823–11832.
9. G. Herrmann, E. Jayamani, G. Mai and W. Buckel, J. Bacteriol., 2008, 190, 784–791.
10. F. Li, J. Hinderberger, H. Seedorf, J. Zhang, W. Buckel and R. K. Thauer, J. Bacteriol., 2008, 190, 843–850.
11. N. P. Chowdhury, A. M. Mowafy, J. K. Demmer, V. Upadhyay, S. Koelzer, E. Jayamani, J. Kahnt, M. Hornung, U. Demmer and U. Ermler, et al. , J. Biol. Chem., 2014, 289, 5145–5157.
12. A.-K. Kaster, J. Moll, K. Parey and R. K. Thauer, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 2981–2986.
13. J. W. Peters, D. N. Beratan, B. Bothner, R. B. Dyer, C. S. Harwood, Z. M. Heiden, R. Hille, A. K. Jones, P. W. King and Y. Lu, et al. , Curr. Opin. Chem. Biol., 2018, 47, 32–38.
14. J. K. Demmer, N. Pal Chowdhury, T. Selmer, U. Ermler and W. Buckel, Nat. Commun., 2017, 8, 1577.
15. J. L. Yuly, P. Zhang, C. E. Lubner, J. W. Peters and D. N. Beratan, Proc. Natl. Acad. Sci. U. S. A., 2020, 117, 21045–21051.
16. C. E. Lubner, D. P. Jennings, D. W. Mulder, G. J. Schut, O. A. Zadvornyy, J. P. Hoben, M. Tokmina-Lukaszewska, L. Berry, D. M. Nguyen and G. L. Lipscomb, et al. , Nat. Chem. Biol., 2017, 13, 655–659.
17. B. G. Galuzzi, A. Mirarchi, E. L. Viganò, L. De Gioia, C. Damiani and F. Arrigoni, J. Chem. Inf. Model., 2022, 62, 4748–4759.
18. J. L. Yuly, P. Zhang, X. Ru, K. Terai, N. Singh and D. N. Beratan, Chem, 2021, 7, 1870–1886.
19. D. H. Evans, Chem. Rev., 2008, 108, 2113–2144.
20. M. C. McCormick, K. Keijzer, A. Polavarapu, F. A. Schultz and M.-H. Baik, J. Am. Chem. Soc., 2014, 136, 8992–9000.
21. R. L. Lord, C. K. Schauer, F. A. Schultz and M.-H. Baik, J. Am. Chem. Soc., 2011, 133, 18234–18242.
22. M.-H. Baik, C. K. Schauer and T. Ziegler, J. Am. Chem. Soc., 2002, 124, 11167–11181.
23. D. H. Evans and K. Hu, J. Chem. Soc., Faraday Trans., 1996, 92, 3983–3990.
24. J. Hu, W. Chuenchor and S. E. Rokita, J. Biol. Chem., 2015, 290, 590–600.
25. Z. Zhou and R. P. Swenson, Biochemistry, 1995, 34, 3183–3192.
26. S. Alagaratnam, G. Van Pouderoyen, T. Pijning, B. W. Dijkstra, D. Cavazzini, G. L. Rossi, W. M. Van Dongen, C. P. Van Mierlo, W. J. Van Berkel and G. W. Canters, Protein Sci., 2005, 14, 2284–2295.
27. D. M. Hoover, C. L. Drennan, A. L. Metzger, C. Osborne, C. H. Weber, K. A. Pattridge and M. L. Ludwig, J. Mol. Biol., 1999, 294, 725–743.
28. P. A. O'Farrell, M. A. Walsh, A. A. McCarthy, T. M. Higgins, G. Voordouw and S. G. Mayhew, Biochemistry, 1998, 37, 8405–8416.
29. Z. Zhou and R. P. Swenson, Biochemistry, 1996, 35, 15980–15988.
30. M. L. Ludwig, K. A. Pattridge, A. L. Metzger, M. M. Dixon, M. Eren, Y. Feng and R. P. Swenson, Biochemistry, 1997, 36, 1259–1280.
31. R. P. Simondsen and G. Tollin, Mol. Cell. Biochem., 1980, 33, 13–24.
32. W. W. Smith, R. M. Burnett, G. D. Darling and M. L. Ludwig, J. Mol. Biol., 1977, 117, 195–225.
33. S. G. Mayhew and M. L. Ludwig, The Enzymes, Elsevier, 1975, vol. 12, pp. 57–118.
34. J. Vervoort, D. Heering, S. Peelen and W. Van Berkel, Methods Enzymol., 1994, 243, 188–03.
35. J. L. Yuly, PhD thesis, Duke University, 2021.
36. N. Mohamed-Raseek and A.-F. Miller, J. Biol. Chem., 2022, 298, 101733.
37. J. K. Demmer, H. Huang, S. Wang, U. Demmer, R. K. Thauer and U. Ermler, J. Biol. Chem., 2015, 290, 21985–21995.
38. Y. Song, J. Mao and M. R. Gunner, J. Comput. Chem., 2009, 30, 2231–2247.
39. D. Bashford and M. Karplus, Biochemistry, 1990, 29, 10219–10225.
40. E. Brunk and U. Rothlisberger, Chem. Rev., 2015, 115, 6217–6263.
41. H. Lin and D. G. Truhlar, Theor. Chem. Acc., 2007, 117, 185–199.
42. H. M. Senn and W. Thiel, Angew. Chem., Int. Ed., 2009, 48, 1198–1229.
43. J. Cheng, X. Liu, J. VandeVondele, M. Sulpizi and M. Sprik, Acc. Chem. Res., 2014, 47, 3522–3529.
44. A. Warshel and M. Karplus, J. Am. Chem. Soc., 1972, 94, 5612–5625.
45. G. Li, X. Zhang and Q. Cui, J. Phys. Chem. B, 2003, 107, 8643–8653.
46. J. P. Götze and M. Bühl, J. Phys. Chem. B, 2016, 120, 9265–9276.
47. H. Jin, P. Goyal, A. K. Das, M. Gaus, M. Meuwly and Q. Cui, J. Phys. Chem. B, 2016, 120, 1894–1910.
48. Y. Zhang, H. Liu and W. Yang, J. Chem. Phys., 2000, 112, 3483–3492.
49. T. H. Rod and U. Ryde, Phys. Rev. Lett., 2005, 94, 138302.
50. O. Carrillo-Parramon, S. Del Galdo, M. Aschi, G. Mancini, A. Amadei and V. Barone, J. Theor. Comput. Chem., 2017, 13, 5506–5514.
51. J. Heimdal, M. Kaukonen, M. Srnec, L. Rulíšek and U. Ryde, Chemphyschem, 2011, 12, 3337–3347.
52. J. Bentzien, R. P. Muller, J. Florián and A. Warshel, J. Phys. Chem. B, 1998, 102, 2293–2301.
53. M. H. Olsson, G. Hong and A. Warshel, J. Am. Chem. Soc., 2003, 125, 5025–5039.
54. H. Li, S. P. Webb, J. Ivanic and J. H. Jensen, J. Am. Chem. Soc., 2004, 126, 8010–8019.
55. M. Cascella, A. Magistrato, I. Tavernelli, P. Carloni and U. Rothlisberger, Proc. Natl. Acad. Sci. U. S. A., 2006, 103, 19641–19646.
56. M. van den Bosch, M. Swart, J. G. Snijders, H. J. Berendsen, A. E. Mark, C. Oostenbrink, W. F. van Gunsteren and G. W. Canters, Chembiochem, 2005, 6, 738–746.
57. L. Shen, X. Zeng, H. Hu, X. Hu and W. Yang, J. Theor. Comput. Chem., 2018, 14, 4948–4957.
58. M. Sundararajan, I. H. Hillier and N. A. Burton, J. Phys. Chem. A, 2006, 110, 785–790.
59. R.-Z. Liao, J.-X. Zhang, Z. Lin and P. E. Siegbahn, J. Catal., 2021, 398, 67–75.
60. C. Greco, M. Bruschi, P. Fantucci, U. Ryde and L. De Gioia, J. Am. Chem. Soc., 2011, 133, 18742–18749.
61. R. N. Tazhigulov, P. K. Gurunathan, Y. Kim, L. V. Slipchenko and K. B. Bravaya, Phys. Chem. Chem. Phys., 2019, 21, 11642–11650.
62. P. Saura and V. R. Kaila, J. Am. Chem. Soc., 2019, 141, 5710–5719.
63. S. Bhattacharyya, M. T. Stankovich, D. G. Truhlar and J. Gao, J. Phys. Chem. A, 2007, 111, 5729–5742.
64. C. Cheng and S. Hayashi, J. Theor. Comput. Chem., 2021, 17, 1194–1207.
65. A. Churg and A. Warshel, Biochemistry, 1986, 25, 1675–1681.
66. M. Breuer, P. Zarzycki, J. Blumberger and K. M. Rosso, J. Am. Chem. Soc., 2012, 134, 9868–9871.
67. I. Daidone, L. Paltrinieri, A. Amadei, G. Battistuzzi, M. Sola, M. Borsari and C. A. Bortolotti, J. Phys. Chem. B, 2014, 118, 7554–7560.
68. J. C. Schöneboom, H. Lin, N. Reuter, W. Thiel, S. Cohen, F. Ogliaro and S. Shaik, J. Am. Chem. Soc., 2002, 124, 8142–8151.
69. J. C. Schöneboom, F. Neese and W. Thiel, J. Am. Chem. Soc., 2005, 127, 5840–5853.
70. C. M. Bathelt, J. Zurek, A. J. Mulholland and J. N. Harvey, J. Am. Chem. Soc., 2005, 127, 12900–12908.
71. S. Shaik, D. Kumar, S. P. de Visser, A. Altun and W. Thiel, Chem. Rev., 2005, 105, 2279–2328.
72. J. N. Harvey, C. M. Bathelt and A. J. Mulholland, J. Comput. Chem., 2006, 27, 1352–1362.
73. J. Mao, K. Hauser and M. Gunner, Biochemistry, 2003, 42, 9829–9840.
74. J.-D. Maréchal, G. Barea, F. Maseras, A. Lledós, L. Mouawad and D. Perahia, J. Comput. Chem., 2000, 21, 282–294.
75. D. E. Bikiel, L. Boechi, L. Capece, A. Crespo, P. M. De Biase, S. Di Lella, M. C. G. Lebrero, M. A. Martí, A. D. Nadra and L. L. Perissinotti, et al. , Phys. Chem. Chem. Phys., 2006, 8, 5611–5628.
76. A. Nemukhin, B. Grigorenko, I. Topol and S. Burt, Int. J. Quantum Chem., 2006, 106, 2184–2190.
77. A. R. Groenhof, M. Swart, A. W. Ehlers and K. Lammertsma, J. Phys. Chem. A, 2005, 109, 3411–3417.
78. P. E. Siegbahn, Biochim. Biophys. Acta, Bioenerg., 2013, 1827, 1003–1019.
79. P. E. Siegbahn, Chem.–Eur. J., 2006, 12, 9217–9227.
80. L. Zanetti-Polzi, C. A. Bortolotti, I. Daidone, M. Aschi, A. Amadei and S. Corni, Org. Biomol. Chem., 2015, 13, 11003–11013.
81. N. J. Fowler, C. F. Blanford, J. Warwicker and S. P. de Visser, Chem.–Eur. J., 2017, 23, 15436–15445.
82. M. Sulpizi, S. Raugei, J. VandeVondele, P. Carloni and M. Sprik, J. Phys. Chem. B, 2007, 111, 3969–3976.
83. L. Castro, L. E. Crawford, A. Mutengwa, J. P. Götze and M. Bühl, Org. Biomol. Chem., 2016, 14, 2385–2389.
84. J. Blumberger, Chem. Rev., 2015, 115, 11191–11238.
85. C. E. Wise, A. E. Ledinina, D. W. Mulder, K. J. Chou, J. W. Peters, P. W. King and C. E. Lubner, Proc. Natl. Acad. Sci. U. S. A., 2022, 119, e2117882119.
86. P. Mondal, K. Schwinn and M. Huix-Rotllant, J. Photochem. Photobiol., A, 2020, 387, 112164.
87. R. K. Kar, A.-F. Miller and M.-A. Mroginski, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2022, 12, e1541.
88. R. E. Georgescu, E. G. Alexov and M. R. Gunner, Biophys. J., 2002, 83, 1731–1748.
89. E. Alexov and M. Gunner, Biophys. J., 1997, 72, 2075–2093.
90. W. D. Cornell, P. Cieplak, C. I. Bayly, I. R. Gould, K. M. Merz, D. M. Ferguson, D. C. Spellmeyer, T. Fox, J. W. Caldwell and P. A. Kollman, J. Am. Chem. Soc., 1995, 117, 5179–5197.
91. C. I. Bayly, P. Cieplak, W. Cornell and P. A. Kollman, J. Phys. Chem., 1993, 97, 10269–10280.
92. M. e. Frisch, G. Trucks, H. B. Schlegel, G. Scuseria, M. Robb, J. Cheeseman, G. Scalmani, V. Barone, G. Petersson, H. Nakatsuji, et al., Gaussian 16, Revision A.03, 2016.
93. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B, 1988, 37, 785–789.
94. M. J. Frisch, J. A. Pople and J. S. Binkley, J. Chem. Phys., 1984, 80, 3265–3269.
95. L. Li, C. Li, Z. Zhang and E. Alexov, J. Chem. Theory Comput., 2013, 9, 2126–2136.
96. L. Wang, L. Li and E. Alexov, Proteins, 2015, 83, 2186–2197.
97. L. Wang, M. Zhang and E. Alexov, Bioinformatics, 2016, 32, 614–615.
98. S. Pahari, L. Sun, S. Basu and E. Alexov, Proteins, 2018, 86, 1277–1283.
99. Y. Song, J. Mao and M. R. Gunner, Biochemistry, 2003, 42, 9875–9888.
100. Y. Song, J. Mao and M. Gunner, Biochemistry, 2006, 45, 7949–7958.
101. The PyMOL MOlecular Graphics System, Version 1.8, Schrödinger, LLC.
102. W. L. DeLano, CCP4 Newsl. Protein Crystallogr., 2002, 40, 82–92.
103. H. M. Berman, J. Westbrook, Z. Feng, G. Gilliland, T. N. Bhat, H. Weissig, I. N. Shindyalov and P. E. Bourne, Nucleic Acids Res., 2000, 28, 235–242.
104. D. L. Nelson, A. L. Lehninger and M. M. Cox, Lehninger principles of biochemistry, Macmillan, 2008.
105. M. J. Harms, C. A. Castañeda, J. L. Schlessman, G. R. Sue, D. G. Isom and B. R. Cannon, et al. , J. Mol. Biol., 2009, 389, 34–47.
106. K. Bartik, C. Redfield and C. M. Dobson, Biophys. J., 1994, 66, 1180–1184.
107. S. Kuramitsu and K. Hamaguchi, J. Biochem., 1980, 87, 1215–1219.
108. T. Takahashi, H. Nakamura and A. Wada, Biopolymers, 1992, 32, 897–909.
109. Y. Y. Sham, I. Muegge and A. Warshel, Proteins, 1999, 36, 484–500.
110. T. J. You and D. Bashford, Biophys. J., 1995, 69, 1721–1733.
111. J. Antosiewicz, J. A. McCammon and M. K. Gilson, Biochemistry, 1996, 35, 7819–7833.
112. C. A. Fitch, G. Platzer, M. Okon, B. Garcia-Moreno E and L. P. McIntosh, Protein Sci., 2015, 24, 752–761.
113. M. J. Harms, J. L. Schlessman, G. R. Sue and B. García-Moreno, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 18954–18959.

---