Open Access Article
This Open Access Article is licensed under a
Creative Commons Attribution 3.0 Unported Licence

Resolving the structure of the E1 state of Mo nitrogenase through Mo and Fe K-edge EXAFS and QM/MM calculations

Casey Van Stappen, Albert Thor Thorhallsson, Laure Decamps, Ragnar Bjornsson* and Serena DeBeer*
Max-Planck Institute for Chemical Energy Conversion, Stiftstrasse 34-36, 45470 Mülheim an der Ruhr, NRW, Germany. E-mail:;

Received 6th May 2019 , Accepted 3rd September 2019

First published on 4th September 2019

Biological nitrogen fixation is predominately accomplished through Mo nitrogenase, which utilizes a complex MoFe7S9C catalytic cluster to reduce N2 to NH3. This cluster requires the accumulation of three to four reducing equivalents prior to binding N2; however, despite decades of research, the intermediate states formed prior to N2 binding are still poorly understood. Herein, we use Mo and Fe K-edge X-ray absorption spectroscopy and QM/MM calculations to investigate the nature of the E1 state, which is formed following the addition of the first reducing equivalent to Mo nitrogenase. By analyzing the extended X-ray absorption fine structure (EXAFS) region, we provide structural insight into the changes that occur in the metal clusters of the protein when forming the E1 state, and use these metrics to assess a variety of possible models of the E1 state. The combination of our experimental and theoretical results supports that formation of E1 involves an Fe-centered reduction combined with the protonation of a belt-sulfide of the cluster. Hence, these results provide critical experiment and computational insight into the mechanism of this important enzyme.


Mo nitrogenase (Mo N2ase) performs a crucial step in the biogeochemical nitrogen cycle, reducing N2 to two molecules of NH3. This enzyme utilizes a two-component system comprised of the active-site containing MoFe protein and reducing Fe protein (FeP). The subunits of the dimeric FeP are connected through a 4Fe–4S cluster, which serves to transfer electrons to the MoFe protein. The MoFe protein is composed of an α2β2 heterotetramer, where each αβ subunit contains two large iron–sulfur clusters, namely the 8Fe–7S P-cluster, which serves as an electron transfer site, and the MoFe7S9C cluster, commonly referred to as the FeMo-cofactor (FeMoco), which serves as the catalytic active site for N2 reduction.1,2

During native turnover, electrons are transferred to FeMoco in a discreet, step-wise fashion. This is initiated by the binding of reduced ATP-bound Fe protein to MoFe, which induces a conformationally gated one-electron transfer from the P-cluster to FeMoco (forming P1+) that is in turn followed by a rapid one-electron transfer from FePred to the P-cluster in a “deficit spending” electron transfer process.3,4 This is followed by hydrolysis of ATP to ADP, the release of two phosphate ions (Pi), and subsequent dissociation of oxidized FeP.5 This process is repeated a total of four times to initiate binding of N2 to the FeMoco cluster and an additional four times to subsequently reduce N2 to 2NH3. In the absence of N2 or other possible substrates, no more than four e/H+ equivalents have been observed to accumulate, and release of H2 leads to relaxation of the cluster to its resting state (Scheme 1).

image file: c9sc02187f-s1.tif
Scheme 1 Abbreviated version of the Lowe–Thorneley cycle emphasizing the states formed in the absence of N2 substrate.6–9

This cycle is commonly interpreted in terms of the pioneering work of Lowe and Thorneley, who described the kinetic relationships between the proposed catalytic intermediates En, in which n represents the number of electrons delivered to the FeMoco cofactors from FeP.6–9 Since each electron transfer step is coupled with the transfer of a proton, H2 may be produced as a side product from any state between n = 2 and 4. To avoid producing excessive amounts of H2, the rate of electron transfer must be fast relative to the rate of H2 production. The electron transfer rate, and in turn population of the various states En prior to N2 binding at the E4 state, may be controlled by varying the ratio of the two protein components. Based on this scheme, it is possible to selectively populate only the E0 and E1 states under conditions in which the rate of H2 production from the E2 state (reported as up to ∼250 s−1)6 is faster than the rate of E2 formation. This is enabled by use of a large ratio of MoFe[thin space (1/6-em)]:[thin space (1/6-em)]FeP (for example, 25[thin space (1/6-em)]:[thin space (1/6-em)]1 or greater) which results in a low rate of electron-transfer. Meanwhile, in the absence of N2, sufficiently low ratios of MoFe[thin space (1/6-em)]:[thin space (1/6-em)]FeP (less than 1[thin space (1/6-em)]:[thin space (1/6-em)]4) should result in the near complete population of E4.10

Recently, significant advances have been made in our understanding of how FeMoco may be capable of accumulating and storing protons and electrons in E2 and E4 as a result of 1H and 57Fe ENDOR studies. These studies support the formation of hydride species in both E2 and E4, which may serve to level the redox potential of the cluster.11–13 However, the E1 and E3 states have remained uncharacterized by EPR methods due in part to their integer spin. In addition, the inability to isolate pure intermediates during turnover has limited the application of other spectroscopic techniques. To-date, there are no reports characterizing the electronic or geometric structure of E3, and only three which have investigated E1.14–16

One of these three investigations of E1 utilized Mo and Fe K-edge EXAFS to structurally characterize this state. Interestingly, significant contractions in the Mo–Fe and Mo–O/N distances (−0.06 Å and −0.07 Å, respectively) were reported relative to the E0 state. Similarly, a contraction of ∼0.05 Å was found in the average short Fe–Fe distances.14 Meanwhile, our recent studies have revealed that formation of E1 involves an Fe-centered reduction, and that Mo remains redox innocent.15 These results appear confounding, as such large contractions in bond distances are non-intuitive for a system undergoing either reduction (in the case of Fe) or no effective oxidation state change (in the case of Mo). In addition, a more recent EXAFS study examining the structural changes of FeMoco upon binding of CO during native turnover (using similar electron-flux conditions as the previous EXAFS study) did not report any bond contractions at Mo when moving from the resting to CO bound state; instead, only minor elongations of the Mo–O/N and Mo–Fe distances were found.17

What is the precise nature of the E1 state? On the basis of previous Mo and Fe K-edge XAS and 57Fe Mössbauer studies, it is clear that Fe is reduced in E1.15,16 However, whether E1 additionally involves protonation of a sulfide or the formation of an iron-hydride (either end-on or bridging) at the FeMoco cluster remains unclear. The permutation space for possible species is quite high, as any of the 7 irons or 9 sulfur sites could hypothetically be protonated. Previous computational investigations have found the bridging sulfides (often referred to as S2B, S3A, S5A) to be the most basic,18,19 while S3B has also been suggested as an initial sulfur protonation site (see Fig. 1 for labelling).19–23 The situation is further complicated by whether the Mo-bound homocitrate is protonated; while computational studies have suggested a protonated hydroxyl group of homocitrate in E0 based on crystal structure comparison,24,25 this information is not currently available for the E1 state.

image file: c9sc02187f-f1.tif
Fig. 1 Labelling scheme used presently to describe the iron and sulfur sites of FeMoco and their relative orientation to homocitrate (red) and residues His442, His195. The coordinating O of homocitrate are labeled as O1 for the carboxylic oxygen and O2 for the hydroxyl group.

To shed light on the nature of the E1 state, we have reinvestigated the Fe and Mo EXAFS of Mo N2ase to elucidate the structural changes which occur at Fe and Mo. The revised experimental changes in this state are further used as a metric for QM/MM DFT calculations to provide further insight into the possible nature of this intermediate, particularly regarding sulfide-protonation vs. iron-hydride formation.

Experimental details

Materials and protein purifications

All reagents were obtained from Sigma-Aldrich (St. Louis, MO) or Fisher Scientific (Fair Lawn, NJ) and were used without further purification. Argon and dinitrogen gases were purchased from Westfalen and passed through an activated copper catalyst to remove any traces of dioxygen before use. Mo nitrogenase of Azotobacter vinelandii was produced as described previously.26 Protein concentrations were determined by Lowry assay.27 The purities of these proteins were >95% based on SDS-PAGE analysis with Coomassie staining. All manipulation of proteins and buffers were performed under an Ar atmosphere.

Preparation of freeze-quenched nitrogenase XAS samples

All XAS samples were prepared under an Ar atmosphere and contained final concentrations of 300 μM MoFe, 12 μM FeP, 50 mM Tris, 200 mM NaCl, 2.5 mM MgCl2, and 10 mM sodium dithionite at pH = 7.5. Turnover state samples were prepared by the addition of an “activating” buffer solution containing 60 mM creatine phosphate, 30 mM MgATP, and 50 units per ml creatine phosphokinase (50% of total final volume) to a solution containing 600 μM MoFe and 24 μM FeP, which was then freeze quenched in liquid N2 after being allowed to react for 120 s. Turnover samples contained final concentrations of 30 mM creatine phosphate, 15 mM MgATP, and 25 units per ml creatine kinase. Although FeP is present in both resting and turnover samples, the Fe present from FeP only accounts for approximately 0.5% of the total Fe in these samples as a 25[thin space (1/6-em)]:[thin space (1/6-em)]1 [MoFe][thin space (1/6-em)]:[thin space (1/6-em)][FeP] ratio was employed.

EPR measurements

EPR measurements were performed to quantify the reduction of resting MoFe, and are detailed in Section S1 of the ESI. An average of 50% reduction in the S = 3/2 E0 signal was found for turnover samples, based upon both relative intensity of the g1 feature at g = 4.3 and spin-integration area (Fig. S1-1). No S = 1/2 or 5/2 signals associated with the one electron oxidized P1+ state of the P-cluster were observed.

X-ray spectroscopic measurements

X-ray absorption measurements of intact nitrogenase MoFe and FeP were obtained at the 9–3 beamline of the Stanford Synchrotron Radiation Lightsource (SSRL). The SPEAR storage ring operated at 3.0 GeV in top-off mode with a ∼500 mA ring current. A liquid N2 cooled double-crystal monochromator with Si(220) crystals at ϕ = 0° was used to select the incoming X-ray energy with an intrinsic resolution (ΔE/E) of ∼0.6 × 10−4, and a Rh-coated mirror was used for harmonic rejection. The X-ray beam size was 1 × 4 mm2 (V × H) at the sample position. A liquid helium flow cryostat was used to maintain a ∼20 K sample environment in order to minimize radiation damage and maintain an inert sample environment. Fluorescence measurements were recorded using a Canberra 100-element Ge monolith solid-state detector. Prior to measurements, each sample was checked for signs of radiation damage by performing subsequent five minute scans over the same sample spot. These tests showed the MoFe protein was stable under X-ray irradiation at the Mo K-edge for >90 minutes, and >70 minutes at the Fe K-edge.

For Mo XAS measurements, the energy of the incoming X-rays was calibrated by simultaneous measurement of a Mo foil and assigning the energy of the maximum of the white line to 20[thin space (1/6-em)]016.4 eV. Full XAS scans were collected by scanning the incident energy from 19[thin space (1/6-em)]780 to 21[thin space (1/6-em)]142 eV. All Fe XAS scans were collected by scanning the incident energy from 6882 to 8093 eV, and calibrated by simultaneous measurement of an Fe foil, with the first inflection point set to 7111.2 eV.

PFY-XAS data processing & statistical analysis

In all experiments, individual scans were normalized to the incident photon flux and averaged using the program Athena from the software package Demeter.28 Further processing of spectra including background subtraction and normalization was also performed using Athena,28 following standard protocols for X-ray spectroscopy described below. Statistical analysis of XAS measurements was performed by normalization of individual scans based on edge area, followed by a calculation of the standard deviation (eqn (1)),
image file: c9sc02187f-t1.tif(1)
where σ is the standard deviation, xi is an individual scan, xav is the average over all scans, and j is the total number of scans. Errors provided for difference spectra were propagated using eqn (2),
image file: c9sc02187f-t2.tif(2)
where σxaa–xbb is the standard deviation of the renormalized spectrum generated by subtraction of fraction xb of spectrum “b” from spectrum “a”. In all cases, xa = 1. Where difference spectra are presented, in which xb = 1, eqn (2) simplifies to:
image file: c9sc02187f-t3.tif(3)

EXAFS fitting was performed using the program Artemis, also of the software package Demeter.28 Possible scattering paths for the EXAFS models were initially determined using FEFF 7.0 in combination with a recent high-resolution crystal structure (PDB ID 3U7Q).29 The structural parameters R (bond distance) and σ2 (bond variance) were allowed to vary during fitting refinement for all measured data. A value of S02 of 1 was used in fitting the Mo EXAFS, and 0.9 in fitting the Fe EXAFS. A single ΔE0 parameter was assigned to all scattering paths at a given edge, and allowed to vary in the refinement of the resting state Mo and Fe K-edge EXAFS. This refined value of ΔE0 was subsequently fixed at these best fit values for further analysis of the turnover and E1 EXAFS data. The Mo and Fe spectra of E1 were generated by multiplying the normalized E0 spectrum by 0.50 (quantified by EPR based on the 50% reduction of the E0 S = 3/2 signal during turnover, Fig. S1-1), and subtracting it from the normalized turnover (E0 + E1) spectrum; the result was renormalized by multiplying by two.

Root-mean-square deviations of the changes in distances of the calculations relative to those determined experimentally were also determined. To do so, we first define the change in distance, ΔR, of path i:

image file: c9sc02187f-t4.tif(4)

We can further define the deviation of the calculated change in distance from that determined experimentally as:

ΔRcalc–expi = (ΔRcalci − ΔRexpi) (5)
where ΔRcalci is the calculated change in a given scattering path, and ΔRexpi is the experimentally observed change in a giving scattering distance, E1–E0. The root mean-square deviation can then be calculated as the normalized sum-of-squares.
image file: c9sc02187f-t5.tif(6)
where j is the total number of paths. To account for the uncertainty in the experimentally determined distances, the RMSD can be weighted using a normal distribution:
image file: c9sc02187f-t6.tif(7)
where χ is the calculated mean, μ is the experimental mean, and σ is the experimental standard deviation.

Mo data processing

Background subtraction and normalization of the averaged Mo EXAFS spectrum was performed using a linear regression for the pre-edge region of 19[thin space (1/6-em)]910–19[thin space (1/6-em)]947 eV, and a quadratic polynomial regression for the post-edge region of 20[thin space (1/6-em)]157–20[thin space (1/6-em)]807 eV. Data were splined from k = 0–17.2 Å−1 using an R-background of 1.0 Å and k-weight of 2. The resulting spectrum was k3-weighted to emphasize the high importance of data at higher k. A k-range of 2–16.5 Å−1 was used for the curve fitting analysis of E0 and E0 + E1 giving a maximum resolution of ΔR = 0.108 Å. A k-range of 2–16 Å−1 was used in the curve fitting analysis of the E1 spectrum, giving a maximum resolution of ΔR = 0.112 Å. All data were fit in R-space using an R-range of 1.5 to 3.5 Å. Due to the already considerable complexity of the EXAFS of the MoFe protein, fitting was limited to include only single scattering paths. No smoothing was used at any point in any of the data processing.

Fe data processing

The Fe EXAFS was processed in a similar fashion to that of the Mo EXAFS. Background subtraction and normalization was performed using a linear regression for the pre-edge region of 6990–7005 eV, and a quadratic polynomial regression for the post-edge region of 7160–8200 eV. Data were splined from k = 0–15.9 Å−1 using an R-background of 1.0 and k-weight of 2. A k-range of 2–15.8 Å−1 was used in the curve fitting analysis of E0 and E0 + E1 to provide a maximum resolution of ΔR = 0.114 Å. Meanwhile, a reduced k-range of 2–13 Å−1 was used for E1, with a maximum ΔR = 0.143 Å. All data were fit in R-space using an R-range of 1.5 to 4.0 Å. Due to the already considerable complexity of the EXAFS of the MoFe protein, fitting was limited to include only single scattering paths. In the case of the split short Fe–Fe scattering path model, the Debye–Waller factors of the two short Fe–Fe scatterers are fixed to be equivalent to one another as the scatterers are the same identity and in a similar environment. This was done to minimize the number of free parameters. No smoothing was used at any point in any of the data processing.


The QM/MM models for E1 were based on our previous model for the E0 resting state.24 It is a spherical QM/MM model (42 Å radius) centered on the carbide of FeMoco. In the QM/MM geometry optimizations, the active region consists of 1000 atoms and a QM region of 133 atoms. All QM/MM calculations were performed in Chemshell version 3.7 (ref. 30) using the built-in MM code DL_POLY31 with the CHARMM36 (ref. 32 and 33) forcefield and ORCA version 4.0 (ref. 34) as QM code. The QM region contains the FeMoco cofactor, singly protonated homocitrate (unless otherwise mentioned) and the sidechains surrounding the cluster which are believed to be critical to describing the coordination, asymmetry, and hydrogen-bonding environment around FeMoco. This includes residues directly coordinating FeMoco (His442, Cys275), neighboring charged residues (Arg96, Arg359), those capable of participating in hydrogen bonding (His195, Gln191, Ser278, Glu380), as well as spatially close residues (Val70, Phe381). For further details, see Fig. S6.4-2 of the ESI. All QM/MM calculations used electrostatic embedding, and the link atom scheme with charge-shifting was used to terminate the QM–MM border as implemented in Chemshell.30 The QM calculations used the TPSSh35 hybrid density functional (previously found to describe the cofactor well24), a ZORA scalar relativistic Hamiltonian,36,37 the relativistically recontracted def2-TZVP38–40 basis set on all metal and sulfur atoms, as well as on the homocitrate, carbide, and two H atoms (def2-SVP on other atoms) and the D3 dispersion correction with Becke–Johnson damping.41–43 The RIJCOSX approximation44–47 was used to speed up computation of Coulomb and HF Exchange integrals. Different electronic states of the cofactor for both the E0 model and the different E1 models were explored by the use of broken-symmetry DFT methodology, as used in previous studies by us.24,48,49 This involves flipping the spin on different Fe atoms starting from a high-spin MS = 35/2 determinant (for E0) or MS = 34/2 (for E1) and then converging to antiferromagnetic low-spin states with a particular MS value. This results in different electronic states (broken-symmetry states), that we label according to which Fe ions are “spin-down”. The three lowest energy states (for both E0 and E1) correspond to the “BS7” category by Lovell, et al.50 Section 6.1 of the ESI provides more detail on the nature of these different states.



The normalized Mo and Fe K-edge XANES spectra of MoFe in the resting (E0), turnover (E0 + E1), and E1 states are provided in Fig. S5-1 and S5-2. We have previously reported on the Mo and Fe K-edges of MoFe in the E1 state in detail.15 The changes observed presently at both Mo and Fe are consistent with our previous findings. No significant shifts are found at the Mo K-edge moving from E0 to E1, indicating no change in either oxidation state or coordination at Mo. Meanwhile, at the Fe K-edge, reduced pre-edge and edge intensities are followed by a concomitant increase in intensity at the white-line when moving from E0 to E1, indicative of an Fe-centered reduction.15,48


The k3-weighted Mo EXAFS of the resting E0, turnover (E0 + E1), and E1 states of MoFe are shown in Fig. 2 along with their corresponding Fourier transforms (FT), while fit parameters are provided in Table 1.
image file: c9sc02187f-f2.tif
Fig. 2 (Top) k3-Weighted and (bottom) FT-spectra of the Mo EXAFS of resting (E0), turnover (E0 + E1), and E1 state. Solid black lines denote experimental spectra while red dashed lines indicate best fits using the 3 scattering path model. The FT-spectra are phase corrected for the Mo–S scattering path.
Table 1 Summary of Mo EXAFS fitting parameters for the resting (E0), turnover (E0 + E1), and E1 spectra using both 3 scattering path (Mo–O/N, Mo–S, and Mo–Fe) and 4 scattering path (Mo–O, Mo–N, Mo–S, and Mo–Fe) models. Standard errors are provided for σ2 and R in parentheses as determined from the fitting procedure
Mo FEFF fits, E0 = 20[thin space (1/6-em)]013.895 eV
Sample Path = 3 N S02 σ2 (10−3 Å2) R (Å) k−1) R-Factor
Resting O/N 3 1 4.20 (1.01) 2.217 (0.009) 2–16.5 min. ΔR 0.108 0.0061
S 3 2.42 (0.30) 2.362 (0.003)
Fe-short 3 3.25 (0.18) 2.689 (0.002)
Turnover (E0 + E1) O/N 3 1 3.52 (1.01) 2.221 (0.009) 2–16.5 min. ΔR 0.108 0.0065
S 3 2.57 (0.36) 2.361 (0.004)
Fe-short 3 3.24 (0.19) 2.688 (0.002)
Turnover (E1) O/N 3 1 2.34 (1.31) 2.221 (0.010) 2–16 min. ΔR 0.112 0.0117
S 3 2.81 (0.68) 2.365 (0.005)
Fe-short 3 3.16 (0.29) 2.697 (0.003)

Sample Path = 4 N S02 σ2 (10−3 Å2) R (Å) k−1) R-Factor
Resting O 2 1 2.36 (0.67) 2.198 (0.007) 1.6–16.4 min. ΔR 0.106 0.0024
N 1 1.72 (1.09) 2.303 (0.011)
S 3 2.51 (0.20) 2.365 (0.004)
Fe-short 3 3.18 (0.11) 2.693 (0.001)
Turnover (E0 + E1) O 2 1 2.53 (0.68) 2.197 (0.006) 1.6–16.4 min. ΔR 0.106 0.0024
N 1 2.01 (1.35) 2.309 (0.011)
S 3 2.24 (0.16) 2.361 (0.002)
Fe-short 3 3.16 (0.11) 2.688 (0.001)
Turnover (E1) O 2 1 0.93 (0.77) 2.202 (0.011) 1.6–16 min. ΔR 0.109 0.0108
N 1 3.00 (0.02) 2.302 (0.006)
S 3 2.57 (0.52) 2.367 (0.006)
Fe-short 3 3.18 (0.26) 2.697 (0.003)

No significant changes are observed from a simple comparison of the resting and turnover states. A slight broadening in the |FT| of E1 is seen, though this is likely due to the increase in noise which naturally results from the spectra subtraction. All spectra exhibit two clear shells (∼2.3 and 2.7 Å), as well as a small feature ∼5 Å, previously reported as a long-distance Mo–Fe scatterer.17

X-ray diffraction crystallography has clearly shown that the resting (E0) MoFe protein contains a single unique Mo, which is coordinated by homocitrate, histidine, and three inorganic sulfides from the FeMoco cluster.29,51 Additionally, there are three neighboring Fe atoms at ∼2.7 Å, ranging from 2.67 to 2.73 Å. Due to the relatively high symmetry of the Mo–S and Mo–Fe distances in the cluster, there are four possible short-range single scattering pathways from the Mo site, which include Mo–S, Mo–Fe, Mo–O, and Mo–N (Fig. 3). In a recent high-resolution crystal structure, the average Mo–O distance is approximately 2.18 Å (at 2.16 and 2.19 Å), while that of Mo–N is 2.33 Å.29 While N and O are generally indistinguishable by EXAFS due to their similar mass, the relatively large deviation between Mo–O and Mo–N bond distances determined by this structure (ΔR = 0.15 Å) in combination with the relatively high resolution of the present experiment (minimum ΔR = 0.108 Å for E0 and 0.112 Å for E0 + E1) warranted further investigation. From Fig. 4, it is clear that the Mo–O and Mo–N paths can be distinctly fit for the E0 state, meaning their fit paths are neither entirely destructive or constructive. However, the use of separate fitting paths results in only very minor statistical improvement (Table 1, 3- vs. 4-path fits for E0). Due to the small magnitude of the Mo–O/N path(s) contribution to the overall spectra, the fit parameters for these paths naturally have a larger degree of uncertainty than the heavier S and Fe scatterers.

image file: c9sc02187f-f3.tif
Fig. 3 Short single-scattering paths of Mo in FeMoco. Model based on coordinates obtained from the XRD structure, PDB ID: 3U7Q.29 These include the Mo–N (blue), Mo–O (red), Mo–S (yellow), and Mo–Fe (rust) distances.

image file: c9sc02187f-f4.tif
Fig. 4 k3-Weighted fits of the resting (E0) Mo K-edge EXAFS spectrum using (top) 3 paths and (bottom) 4 paths.

Moving to the reduced E0 + E1 and E1 spectra, we find that no significant changes are observed in the Mo–S distance in either fit model relative to the resting E0 (Fig. 5). A small contraction in the average Mo–Fe distance is found when using either model to fit E0 + E1, while a small expansion is found for E1. As any real change in distance should trend when moving from E0 + E1 to E1 (as in an increase for E0 + E1 should be even larger in E1), we conclude that no significant structural changes occur from the perspective of Mo when progressing from the E0 to E1 state of MoFe.

image file: c9sc02187f-f5.tif
Fig. 5 Comparison of variation in determined bond distances from Mo K-edge EXAFS. Changes in distances are calculated by subtraction of the fit E0 distances from those of the turnover (E0 + E1) and E1 fits. Error bars are reported on a 95% confidence interval.


The major single-scattering paths in FeMoco from the perspective of Fe include Fe–S, Fe–Mo, Fe–Fe(short), and Fe–Fe(long). Additionally, it is known that a central carbide exists in the core of the FeMoco cluster, adding the possibility of a Fe–C path (Fig. 6a and b).29,52,53 Based on the XRD structure (PDB ID: 3U7Q), the FeMoco cluster is relatively symmetric and as a result the deviations in distance for a given path are quite small (for example, ∼0.04 Å in the long Fe–Fe path and ∼0.1 Å in the short Fe–Fe).29 However, the P-cluster presents a more complicated situation. While there are only Fe–S and Fe–Fe single scattering paths to consider for this cluster, its greater asymmetry results in a much wider distribution of distances. This becomes clear when examining distances determined in the high resolution crystal structure (PDB ID: 3U7Q),29 where the short Fe–Fe distances range from 2.50 to 2.92 Å (Fig. S2-2). Similarly, the Fe–S distances vary from 2.25 to 2.47 Å, and the long Fe–Fe path from 3.79 to 5.48 Å. As a result of this distribution of distances, several considerations on how to appropriately treat these three paths in our model must be made, while still maintaining the minimal necessary number of variables. The Fe–S shell may still be modeled by a single path; however, the relatively large variation in these distances results in a greater degree of static disorder, and in turn a larger Debye–Waller factor. While the long Fe–Fe distances in the FeMoco cluster are nearly identical (Fig. 7a), fitting the highly disordered long Fe–Fe distances in the P-cluster (Fig. 7b) would necessitate separately treating numerous long Fe–Fe scattering paths. This approach would require the use of a large number of parameters, and would only make minor contributions to the long-range region. Therefore, the contribution of the long Fe–Fe scattering paths arising from the P-cluster are not considered in our model. The disorder of the short Fe–Fe distances in the P-cluster present an intermediate case, where they are disordered but still tightly clustered enough that they must be considered in the model. There is an existing precedence for (a) splitting the Fe–Fe(short) scattering path into two shells or (b) using a reduced number (N) of short Fe–Fe scatterers.14,54,55 The use of a reduced number N in the short Fe–Fe path operates under the assumption that the individual short Fe–Fe scattering paths in the P-cluster effectively cancel one another, as evidenced by an almost absent short Fe–Fe shell in the FT-EXAFS of the P-cluster only variant ΔnifB.54
image file: c9sc02187f-f6.tif
Fig. 6 Depiction of the unique single-scattering paths for each of the three classes of Fe (a–c) in FeMoco. Generated using coordinates of XRD structure PDB ID 3U7Q.29 These include the Fe–Mo (cyan), Fe–C (gray), Fe–S (yellow), Fe–Fe(short) (rust), and Fe–Fe(long) (red) paths.

image file: c9sc02187f-f7.tif
Fig. 7 |FT| space spectrum of the Fe K-edge EXAFS of E0 (black, solid), scaled from 2.8–6 Å to compare the cumulative Fe–Fe(long) scattering paths for (a) FeMoco and (b) P-cluster, as determined by FEFF calculations of the 3U7Q XRD structure.29 The paths for each unique long Fe–Fe scatterer are depicted in blue (dotted), and the sum of these unique long Fe–Fe scatterers is shown in red (dashed). The high symmetry of the FeMoco cluster results in near-identical Fe–Fe(long) distances, which accumulate to provide a significant contribution to the spectrum, which can be modeled by a single Fe–Fe scattering path. Meanwhile, the long Fe–Fe distances in the P-cluster are highly disordered, making little contribution to the overall fit.

To cover both approaches, models using a single Fe–Fe(short) scattering path, split Fe–Fe(short) scattering paths, and a single Fe–Fe(short) scattering path with reduced N = 1.65 were generated (fit parameters and figures are provided in the ESI). In short, both the use of a split Fe–Fe(short) scattering path and single Fe–Fe(short) path with reduced N can be used to reasonably fit the data, while use of full N with only one Fe–Fe(short) scatterer results in a poor fit. However, use of a split Fe–Fe(short) shell demonstrates that the two Fe–Fe paths are not entirely destructive (Fig. S3-1), and produce a nominal statistical improvement of the fit (Tables S5-3 and S5-6). Therefore, the model used to fit the data presented here utilizes Fe–S, Fe–Mo, Fe–C, Fe–Fe(long), and split Fe–Fe(short) single-scattering paths. While previous EXAFS studies of the FeMoco cluster have found statistical improvement by the inclusion of a light Fe–X scatterer,53 our fits show that inclusion of the Fe–C path has little to no impact on either the statistics of the fit or the parameters determined for the other scattering paths. This discrepancy may in part be due to the presence of the P-cluster in the current samples, which significantly reduces the contribution of the Fe–C path to the overall EXAFS. The feasibility of objectively fitting the Fe–C scattering path is further discussed in Section S3 of the ESI. Despite the lack of significant statistical improvement achieved by inclusion of the Fe–C scatterer, we have opted to still include this path in our presented 6-path model in acknowledgement of its presence.

The k3-weighted Fe EXAFS spectra of resting (E0), turnover (E0 + E1), and E1 MoFe are provided in Fig. 8 along with the corresponding FTs and best fits using the 6-component model. The corresponding parameters for these fits are provided in Table 2.

image file: c9sc02187f-f8.tif
Fig. 8 Fe K-edge EXAFS spectra of resting (E0), turnover (E0 + E1), and E1 state in both (top) k-space and (bottom) |FT| space. The k-space spectra are k3-weighted phase-shifted and the |FT| are phase corrected for the Fe–S scattering path. Solid black lines denote experimental spectra while red dashed lines indicate best fits using the 6 scattering path model.
Table 2 Fe K-edge EXAFS fitting parameters of the 6 component model as applied to the resting (E0), turnover (E0 + E1), and E1 spectra. Statistical errors are provided for σ2 and R as determined from the fitting procedure
Fe FEFF fits, E0 = 7118.03 eV
Sample Path = 6 Na S02 σ2 (10−3 Å2) R (Å) k−1) R-Factor
a Please see Section S2 of the ESI for details on the determination of the path degeneracy N.
Resting C 0.4 0.9 2.00 (0.82) 1.911 (0.034) 2–15.8 min. ΔR 0.114 0.0142
Mo 0.2 2.00 (0.62) 2.683 (0.013)
S 3.6 5.59 (0.44) 2.267 (0.002)
Fe-short 2.53 5.72 (1.30) 2.622 (0.003)
Fe-short 0.93 2.854 (0.012)
Fe-long 0.8 1.99 (0.09) 3.691 (0.010)
Turnover (E0 + E1) C 0.4 0.9 2.00 (2.81) 1.928 (0.024) 2–15.8 min. ΔR 0.114 0.0093
Mo 0.2 2.00 (2.27) 2.679 (0.010)
S 3.6 6.10 (0.39) 2.273 (0.001)
Fe-short 2.53 6.45 (1.11) 2.620 (0.003)
Fe-short 0.93 2.858 (0.010)
Fe-long 0.8 1.99 (0.85) 3.704 (0.007)
Turnover (E1) C 0.4 0.9 2.00 (0.35) 1.992 (0.029) 2–13 min. ΔR 0.143 0.0079
Mo 0.2 5.00 (0.00) 2.664 (0.001)
S 3.6 7.12 (0.00) 2.283 (0.003)
Fe-short 2.53 7.60 (0.91) 2.603 (0.003)
Fe-short 0.93 2.856 (0.010)
Fe-long 0.8 4.66 (2.10) 3.713 (0.014)

Interestingly, the splitting between the Fe–S and Fe–Fe(short) shells decreases moving from E0 to E0 + E1, and is completely lost in the simulated pure E1 spectrum. Fitting these data, it is apparent from the Debye–Waller factors of the Fe–S and Fe–Fe(short) shells that the degree of disorder of the metal clusters of MoFe increase when under turnover (Table 2 and Fig. 9). As all spectra were measured at the same temperature, it is sensible that this effect likely arises from static disorder. This is expected, considering the E0 + E1 spectrum represents a 50/50 mixture of the E0 and E1 states.

image file: c9sc02187f-f9.tif
Fig. 9 Comparison of σ2 (left) and |FT| (right) of the Fe–S, Fe–Fe(short 1), and Fe–Fe(short 2) paths for the 6-component model fits of resting (E0), turnover (E0 + E1), and E1 MoFe.

The similarity of the Fe–Mo with Fe–Fe(short) distances combined with the relatively small contribution of the Fe–Mo scattering path to the overall fit results in considerable correlation between the fit Fe–Mo and Fe–Fe(short) parameters when determined from Fe EXAFS. As a result, there is an intrinsically large degree of error in the Fe EXAFS determined Fe–Mo distances. Regardless of the applied model, no significant changes in the Fe–Mo are observed, consistent with results obtained from the Mo EXAFS. Generally, a small but statistically significant increase in the Fe–S distances is found, while a small decrease in the Fe–Fe(short 1) distance is observed in fitting the E1 spectrum (Fig. 10). These changes hold true for all models investigated in the present study (Fig. S5-9). It is worthwhile to note that while the Fe–S path represents an average of both FeMoco and P-cluster, the Fe–Fe(long) path is predominately representative of the FeMoco cluster. Additionally, by splitting the short Fe–Fe path, the first Fe–Fe(short 1) path represents a combination of the FeMoco and P-cluster Fe–Fe(short) scatterers, while Fe–Fe(short 2) should only be representative of the P-cluster. It is therefore not surprising to see that Fe–Fe(short 1) contracts while no effective changes are observed for Fe–Fe(short 2).

image file: c9sc02187f-f10.tif
Fig. 10 Comparison of variation in determined bond distances from the resting E0 state in the 6-component model applied to the Fe EXAFS. Distances are calculated by subtraction of the fit E0 distances from those of the E0 + E1 and E1 fits. The differences in the Fe–Fe(short) distances are split such that Fe–Fe(short 1) is to the left, and Fe–Fe(short 2) to the right. Error bars are reported on a 95% confidence interval.

Thus far we have rigorously reexamined the Mo and Fe EXAFS of the E0 and E1 states of MoFe using a variety of possible models. The presented results have demonstrated that no effective changes are found in the first coordination sphere of the Mo of the FeMoco cluster, while all fit models show small but consistent variations in the Fe distances. With these results in hand, we turn to theoretical methods, specifically QM/MM, for a more detailed investigation of the specific nature of the E1 state.

QM/MM calculations

We performed a series of quantum mechanics/molecular mechanics (QM/MM) calculations to investigate the possible electronic states accessible by single e reduction of FeMoco and the nature and position of the transferred H+ in the FeMoco cluster (i.e. iron-hydride formation vs. sulfur protonation). A previous QM/MM study of resting state MoFe protein by one of us has demonstrated that very good agreement between the computed FeMoco structure and the high-resolution crystal structure can be obtained when the protein environment is included via QM/MM and the TPSSh hybrid functional is used.24 Here, we employ an almost identical computational approach to study the changes which occur going from the E0 to E1 states. To reduce systematic errors present in the calculations (such as any over- or underestimation of the covalency of specific chemical bonds, or errors resulting from the simplified spin-coupling treatment employed), we will focus on the changes which occur in E1 relative to the resting E0 state. Additionally, in order to compare the changes in calculated bond distances ΔR (defined by R(E1) − R(E0) for a given path), with those determined by EXAFS, we have averaged the Mo and Fe bond distances to the same level of resolution as observed in the Mo EXAFS 4-path model, and the Fe EXAFS 6-path model.

In approaching the question of the identity of E1, three major considerations must be made, namely (a) the protonation state of the Mo-binding homocitrate in E0 and E1, (b) the location of protonation on the FeMoco cluster (as well as whether it even occurs), and (c) the possible electronic states accessible upon reduction, directly related to the site of reduction as discussed later. Additionally, as it is unknown whether the Nδ or Nε position of His195 is protonated in E1, we have considered both possibilities.

The protonation state of homocitrate has been previously discussed for the E0 state, where comparison of computed distances (primarily the distance between the oxygens of the hydroxyl and carboxylate groups in homocitrate) and quantum crystallographic refinement have indicated the Mo-bound hydroxyl group is protonated.24,25 We revisit this protonation assignment here in the context of E1, and discuss it in greater detail in Section 6.2 of the ESI. Briefly, all models which leave the homocitrate unprotonated (in either the resting E0 or reduced E1 states) result in the highest occupied molecular orbitals (HOMOs) having positive energies, which is unphysical. Additionally, models involving deprotonated homocitrate result in the contraction of the averaged Mo–O bond by 0.05 Å and lengthening of both Mo–S and Mo–Fe by ∼0.02 Å when moving from E0 to E1, disagreeing with the EXAFS results. Meanwhile, all models which retain a protonated homocitrate in both E0 and E1 states show a mild expansion of the average Mo–O/N distance by ∼0.02 Å, a negligible expansion of the Mo–S distances, and a small contraction of the Mo–Fe distance by ∼−0.02 Å. The combination of the strong Mo–O contraction, which is well out of the standard error of the experiment, combined with the non-physical energies of the HOMOs indicate that the Mo-bound hydroxyl group of homocitrate is likely protonated in both the E0 and E1 states, and will be considered as such for the remainder of the results.

We have found that addition of an electron to the FeMoco in the absence of a proton results in mild modulation of all bond lengths, which are in reasonable agreement with those observed in the EXAFS (Fig. S6.4-3 through S6.4-8). However, similar to the case of the unprotonated homocitrate, the energies of the HOMOs become positive (see Tables S6.4-3 through S6.4-6), suggesting that protonation of FeMoco may be required to obtain a physically relevant E1 state. It has already been established that the belt sulfides of the cluster (S2B, S3A, and S5A) are the most basic according to another detailed QM/MM study,18 and we have therefore focused the present study on the protonation of these sulfurs. Additionally, we have investigated the possibility of terminal hydride formation at the sterically unhindered Fe6, as well as the protonation of S1A and S3B sites. While a bridging hydride model was considered, it was found to be unstable and consistently converted to terminal coordination upon optimization; previous computational studies of several bridging hydride models of E1 were also found to be highly energetically unfavorable.18 The direction of protonation in each of these models is discussed in Section S6.3.2 of the ESI. Finally, we investigated the possibility of a carbide protonation in E1, as this has been suggested to be thermodynamically favorable according to calculations by Rao, et al.56 For each of these protonation states, we have further considered three unique electronic states of the cofactor (referred to as broken-symmetry (BS) solutions BS-235, BS-346 and BS-247) which correspond to different locations of the Fe up/down local spins and delocalized pairs in the cluster. In this nomenclature, the first two numbers indicate the two Fe which form a spin-down mixed-valent delocalized pair in the [Fe4S3C] sub-cubane, and the third denotes the Fe which becomes “spin-down” in the [MoFe3S3C] sub-cubane. This third Fe is spin-localized, is either ferric (in E0) or ferrous (in E1). Section S5.1 of the ESI provides a more detailed description of this scheme. Although all three BS solutions share the same favorable antiferromagnetic spin orientation (known as BS7 in the literature),50 it is necessary to consider all three since the effective C3v symmetry of the cluster is broken by the secondary environment of the protein. In the context of E1, these three solutions effectively rotate which iron in the cluster is reduced (Fe4 or Fe5 in BS-235, Fe2 or Fe6 in BS-346, and Fe3 or Fe7 in BS-247; see Fig. 1). While it is possible that the added electron could be on either of two Fe atoms for a given BS solution, we have found that reduction is always localized to one of the Fe atoms of the [MoFe3S3C] sub-cubane (Fe5 through Fe7). Other electronic states that would hypothetically result in reduction of the [Fe4S3C] cubane were considered but found to be energetically unfavorable (>15 kcal mol−1), in agreement with a previous study by Cao et al.18

The presented broken-symmetry solutions have a final spin-state of MS = 2, which was found to be generally energetically more favorable than MS = 1 (with the Fe6-hydride state being the only exception, where the MS = 1 and MS = 2 energies are comparable). This is the logical spin state when a spin-down high-spin Fe(III) (local spin 5/2) is reduced to spin-down Fe(II) (local spin 2), changing the total spin-state of the cluster from MS = 3/2 (E0) to MS = 2 (E1) within a highly simplified spin-coupling model. This is in good agreement with the experimental spin states of [MoFe3S4]3+ (S = 3/2) and [MoFe3S4]2+ (S = 2) synthetic cubanes.57,58 Further details about the electronic states can be found in the ESI, Section S6.1.

To compare the structures of the E1 models featuring different protonation/electronic states to the new EXAFS data in an unbiased fashion, it is helpful to first introduce and discuss a simple nonbiased metric. In this regard, the root-mean-square deviation (RMSD) of the R(E1) − R(E0) structural differences relative to the experimental EXAFS data is quite useful; however, as the experimental standard deviation for different scattering paths can vary considerably (see Fig. 5 and 10), a Gaussian-based weighting scheme has been employed that takes the experimental deviation into account in the RMSD metric. This approach is detailed in the statistical analysis section of the Experimental details.

Fig. 11 shows the weighted RMSDs for the different computational models, where the results for each protonation state have been averaged over the BS electronic states (for the RMSDs of individual BS solutions, see Fig. S6.4-10 and 6.4-11 in the ESI). The results most clearly reveal that the model involving a protonated carbide strongly deviates from the EXAFS results, featuring an RMSD >0.06 Å. The most favorable models also vary considerably depending on whether the Nε or Nδ positions of His195 are protonated. In the case of His195-Nδ(H), the S2B(H) position appears to most favorably agree, followed by the S1A(H); both the S3A(H) and “no H+” models appear poor. Of the models involving His195-Nε(H), the S5A(H) and S2B(H) are in most favorable agreement, while S3B(H) fits relatively poorly. It is worthwhile to note that for most models the Nε protonated state of His195 generally gives consistently lower RMSDs than the Nδ protonated state.

image file: c9sc02187f-f11.tif
Fig. 11 Weighted root-mean-square deviations for the calculated changes in distances (ΔRcalc) of the computational models of the E1 state averaged across the BS-235, BS-346, and BS-247 solutions. Each individual contributing path has been weighted based on the experimental standard deviation, as described in the statistical analysis section of the Experimental details. Experimental changes in distances (ΔRexp) used in these calculations were acquired from E1–E0.

To distinguish these models further, we turn to how well each reproduces the experimentally observed changes in the individual scattering paths. As an example, the changes in the average calculated Mo and Fe distances is provided in Fig. 12 for the S2B(H) and Fe6(H) models. Section 5.6 of the ESI contains complete comparisons of the calculated changes in bond distances for all E1 models, with protonation of either the Nε or Nδ positions of His195, in all three considered BS solutions, relative to the complimenting three BS solutions of the E0 state. All three BS solutions for all models display an increase in the average Mo–O bond length. This increase in Mo–O bond length is relatively mild in most cases, ranging from 0.01–0.03 Å, although the “no H+” model consistently results in expansions of ∼0.04–0.05 Å. Although the EXAFS cannot resolve the individual Mo–O paths, we note that this calculated increase is almost completely attributable to the Mo–O(2) bond (Fig. 1), which corresponds to the hydroxy group of homocitrate. No effective increase in the Mo–NHis442 distance is found in any of the calculations, for any model or BS solution. All models display a mild increase in the average Mo–S distance; most are within the range of experimental error, although the S3B(H) model presents a significantly greater expansion than the others (∼0.03 Å). Lastly, all belt sulfide models (as well as S1A(H)) show small decreases of 0.005–0.03 Å in the average Mo–Fe distance; this decrease is particularly exacerbated in the case of S3A(H) (Fig. S6.4-5).

image file: c9sc02187f-f12.tif
Fig. 12 Comparison of the change in calculated distances upon reduction from E0 → E1 for the S2B(H), S5A(H), and Fe6(H) models in the BS-235 electronic state, with His195 treated as Nε(H). Experimental Fe–Fe(short) distances correspond with the (left) Fe–Fe(short 1), N = 2.53 and (right) Fe–Fe(short 2), N = 0.93 shells described in Table 2. E0 calculated distances from the BS-235 solution are used in the calculated differences, as previous comparisons have shown this solution to best reproduce those determined in the 3U7Q crystal structure.24 A full comparison of all calculated changes in relevant bond distances for all 6 protonated models (in both the Nε(H) and Nδ(H) singly-protonated states of His195), in all three considered BS-solutions, relative to all three BS solutions of the E0 state are provided in Fig. S5.6-3 through S5.6-8. Error bars for the experimentally determined distances are reported on a 95% confidence interval.

Comparing the average calculated Fe distances with experiment, we find that the expansion of the Fe–S path is reasonably reproduced by the three belt-sulfide protonated models, as well as S1A(H) and “no H+”. Meanwhile, the Fe6-hydride model displays a decrease in average Fe–S distance, counter to our EXAFS results. No model is capable of reproducing the contraction found in the Fe–Fe(short 1) distance, although the three belt-sulfide protonated models show no effective change while all others show an expansion (Fig. S6.4-6 and S6.4-7). The C(H) model is particularly extreme, with an expansion of ∼0.08 Å. All models also predict virtually no change in the average Fe–Fe(long) distance with the exception of C(H), which shows an expansion of up to ∼0.14 Å (Fig. S6.4-8).

From our comparisons with the present EXAFS results, we conclude that the S2B(H), S5A(H), and S1A(H) of models of E1 appear the most reasonable. While the weighted RMSD metric suggested the Fe6(H) hydride model as also reasonable, the considerable overestimation of the change in Mo–O and the wrong trend in the change of average Fe–S distance make this model less likely.

With this experimental parameterization in hand, we can compare the relative energies of these different states (Fig. 13). We find that the S2B(H) model is most favorable when Nδ(H) is used at His195, with the S5A(H) model appearing ∼7 kcal mol−1 higher in energy. However, the energies of the S5A(H) and S2B(H) models become comparable when the Nε position of His195 is protonated. Meanwhile, the S3A(H), S1A(H), Fe6(H), S3B(H), and C(H) models all appear considerably unfavorable in either of the His195 protonation states, with energies >10 kcal mol−1 higher than those of the lowest energy solutions.

image file: c9sc02187f-f13.tif
Fig. 13 Relative energies of the calculated energies of the investigated protonated models of E1 with both (a) Nδ(H) and (b) Nε(H) protonation states of His195. Displayed energies are QM/MM total energies.


The present results show that very little change occurs in the metal clusters of MoFe upon formation of E1, counter to the results of a previous EXAFS study.14 Particularly, we find that the large contractions in Mo–O/N, Mo–Fe, and Fe–Fe(short 1) previously reported from Mo and Fe K-edge EXAFS do not occur in the present system (Table 3). One plausible explanation for the large discrepancies between the present study and the previously reported14 Mo and Fe K-edge EXAFS lies in the fitting procedure used to interpret the data. In particular, the models employed in the previous EXAFS study utilized significantly path-dependent values of ΔE0, the parameter which is used to align the energy grids of the experimental spectrum and the model.14 Fundamentally, E0 (not to be confused with the resting E0 state of N2ase) is used to describe the kinetic energy necessary for a photoelectron to escape an absorber. Therefore, it is characteristic of the absorbing atom and independent of the scattering path. Small variations in ΔE0, on the order of up to a couple eV, are generally acceptable when describing the paths of two unique absorbers of the same type in the same sample. For example, EXAFS fitting of NiFe or FeFe hydrogenases may require two unique ΔE0, one for the Fe of the 4Fe–4S clusters, and a second for the Fe–Fe active site. In the previous report of the Mo and Fe EXAFS of MoFe, large variations in ΔE0 of up to 17 eV and 7 eV were employed for different paths involving Mo and Fe, respectively.14 As only a single species of Mo is present in these samples, there is no physical justification for using radically different values of ΔE0 for the various scattering pathways. These changes in ΔE0 effectively change the phase of the individual fit paths relative to one another, making the bond distances determined by said fit effectively arbitrary (please see Section S4 of the ESI for further discussion). This is an unfortunately common mistake in EXAFS fitting;59–64 to quote Scott Calvin – “There are few mistakes more common in published EXAFS work, and more unambiguously wrong, than publishing fits where every path has a unique E0…”.65
Table 3 Comparison of the present EXAFS determined changes in path lengths (E1–E0) with those previously reported.14 “present” distances were determined from the 3-path Mo and 6-path Fe fits. All standard errors are rounded to the nearest third decimal
Path ΔR(E1 − E0) (Å)
Present Ref. 14
a Mo-Fe path distances shown were determined from Mo K-edge EXAFS.
Mo–O/N 0.004 ± 0.014 −0.07
Mo–S 0.003 ± 0.006 0.00
Mo–Fea 0.008 ± 0.004 −0.06
Fe–S 0.016 ± 0.004 0.02
Fe–Fe(short 1) −0.019 ± 0.005 −0.04 to −0.06
Fe–Fe(short 2) 0.002 ± 0.016 −0.01 to −0.02
Fe–Fe(long) 0.022 ± 0.017 0.01 to 0.02

With the experimental bond distances of E0 and E1 in hand, we can compare the observed differences between these two states with those of model complexes and other FeS cluster proteins. Fig. 14 provides the changes in average Fe–S and Fe–Fe distances between several oxidation states in a variety of model cubane clusters and in FeP, based on various literature reports. On average, the Fe–S distance increases by 0.015–0.03 Å in both the cubane model complexes and FeP. Meanwhile, a spread of changes are observed in the Fe–Fe(short) distances, ranging from expansion to contraction. The model cubane complexes show a small decrease on average, while those observed for the FeP are dependent on the redox couple (where a slight expansion is seen going from 2FeII2FeIII → 3FeIIFeIII, and a contraction going from 3FeIIFeIII → 4FeII). We can conclude from these comparisons that both the sign and magnitude of the observed changes are consistent with an Fe-centered reduction.

image file: c9sc02187f-f14.tif
Fig. 14 Comparison of the changes in observed core Fe–S and Fe–Fe distances upon reduction for the turnover (E0 + E1) and E1 states (black and red squares, respectively) with those observed in a series of cubane model complexes57,66–71 (circles) and in FeP72–75 (triangles) based on XRD. The differences in Fe–Fe(short) distances are split to represent Fe–Fe(short 1) to the outside-left, and Fe–Fe(short 2) to the outside-right.

We can now use our revised understanding of the structural changes that occur upon reduction of MoFe from E0 to E1 as a metric to judge the various QM/MM models of E1. The structural changes which occur in several of the QM/MM models when moving from E0 → E1 provide reasonable agreement with the EXAFS data, while several others do not. In nearly all cases, the Mo–O bond length changes were found to be overestimated compared to the EXAFS data. This overestimation relates to the proton of the hydroxy group, which forms a strong hydrogen bond with the carboxylate-arm of the homocitrate (see Fig. S6.4-2). The strength of this hydrogen bond is highly dependent on how well the hydrogen-bonding environment is described in the model, and hence can directly impact the calculated Mo–O distance. We have previously observed spontaneous proton transfer to the carboxylate group in some FeMoco models, revealing a sensitivity to the precise location of this proton; under experimental conditions this acidic proton may even be delocalized between these oxygen atoms, a feature that our computational models are currently unable to capture.

The energetic comparison in Fig. 13 reveals large differences between the investigated models, which are overall in agreement with the result of a previous computational study.18 The bridging sulfides are generally more basic than other sulfides, with S2B and S5A being much more favorable than S3A. The variation between the three belt-sulfide positions is likely related to the peptide backbone environment near S3A, as peptide NH groups show weak hydrogen bonds surrounding this sulfide which render protonation unfavorable. Calculations of these E1 models in the absence of the MM environment (using instead a polarizable continuum model) confirm the secondary coordination environment as the source of this disparity, where protonation becomes almost equally favorable for all belt-sulfide sites (see Fig. S6.3-1 of the ESI). Additionally, the formation of Fe6(H), a terminal hydride, is not thermodynamically favorable at the redox level of E1.

Several studies have suggested that His195 functions as a competent proton donor to the FeMoco cluster;76–78 indeed, the orientation of His195 and 3.2 Å S2B-Nε distance in the 3U7Q crystal structure of E0 suggests that His195 is likely protonated at Nε, forming a hydrogen bond to S2B.29,79 Proton transfer from Nε(H) of His195 to S2B could plausibly result in immediate reprotonation of His195. One possible proton pathway involving His195, Tyr281 and a water molecule has been previously discussed by Dance, where reprotonation of His195 occurs at the Nδ position.19,76 Although this is a plausible scenario, there is no direct experimental evidence for this inverse protonation state in any En state. If a regular protonation state of His195 is considered instead, the S2B and S5A become equally plausible protonation sites in E1, based on the relative energies (∼2 kcal mol−1, see Fig. 13), likely due to reduced basicity of the S2B site via the His195-S2B hydrogen bond. We note that a previous QM/MM study found a larger energy difference (7.6 kcal mol−1) between the S2B(H) and S5A(H) E1 models; the reasons for this difference are not clear but may be related to the use of different functionals (that has been the subject of discussion in the literature80) or slightly different QM/MM models.

As previously mentioned, the electronic state considered appears to determine the specific site of reduction on FeMoco (Fe5, Fe6 or Fe7). One might imagine the electron and proton ending up at the same or neighboring sites (e.g. a redox event at Fe6 resulting in protonation of Fe6-bound S2B, a E0(BS-346) → E1(346)-S2B(H) scenario). Comparing the QM/MM determined energies, we can see that the BS solution which places the additional electron either directly on or neighboring the site of protonation becomes most favorable in all cases. However, this effect is only a few kcal mol−1 at best. This implies that while the basicity/hydricity of the protonation sites investigated are, to a certain degree, sensitive to the location of the additional electron, the surrounding secondary environment of the cluster plays a dominant role in determining the most favorable protonation site.

To conclude, E1 models featuring a protonated bridging sulfide fit best with the EXAFS data. Based on the relative energies of these various models, those involving a protonated sulfide at either the S2B or S5A positions in tandem with an Fe-centered reduction at the [MoFe3S3C] sub-cubane are the most likely candidates to describe the E1 state. These results are in good agreement with numerous experimental studies of ligand-bound states of both Mo and V nitrogenase in which substitution of the belt sulfides occurs.81–85 Further distinguishing between S2B and S5A as protonation site in E1 likely requires experimentally establishing whether His195 can serve as a direct proton donor in the E0 → E1 process. The similar basicity of these sulfide sites also suggest that both may play a role in the formation of other reduced states of FeMoco, e.g. E4, where potentially 4 protons have been added to the cofactor.


Determining how nitrogenases are capable of accumulating 3–4 reducing equivalents, while maintaining effectively the same reduction potential, is critical in understanding how these enzymes bind and activate N2. The E1 state represents the first critical step in this process. The present study has used Mo and Fe K-edge EXAFS to characterize the structural changes which occur upon reduction of E0 to E1. While a previous report claimed large contractions occur in Mo–O, Mo–Fe, and short Fe–Fe distances,14 we have found that only minor modulation of these distances occurs. By comparing our observations with both known model complexes and FeP, we have found that formation of the E1 state is consistent with an Fe-centered reduction, in agreement with our previous observations.15 Furthermore, the combination of the present EXAFS results with QM/MM calculations supports that protonation of a belt sulfide likely occurs in E1, most favorably at the S2B or S5A positions.

Conflicts of interest

There are no conflicts to declare.


C. V. S., L. D., R. B., and S. D. would like to thank the Max-Planck Society for funding. S. D. acknowledges the European Research Council under the European Union's Seventh Framework Programme (FP/2007-2013) ERC Grant Agreement number 615414, and both S.D. and C.V.S. acknowledge the DFG SPP 1927 “Iron–Sulfur for Life” (project DE 1877/1-1) for funding. Additionally, R. B. would like to thank the Icelandic Research Fund grant no. 162880051 and C. V. S. the IMPRS-RECHARGE for further funding. L. D. would like to thank the Peter und Traudl Engelhorn-Stiftung for funding. Additionally, Alexander Gutt and Natalia Levin are thanked for their assistance in the collection of XAS data. Open Access funding provided by the Max Planck Society.


  1. B. K. Burgess and D. J. Lowe, Chem. Rev., 1996, 96, 2983–3012 CrossRef CAS PubMed.
  2. B. M. Hoffman, D. Lukoyanov, Z.-Y. Yang, D. R. Dean and L. C. Seefeldt, Chem. Rev., 2014, 114, 4041–4062 CrossRef CAS PubMed.
  3. K. Danyal, D. R. Dean, B. M. Hoffman and L. C. Seefeldt, Biochemistry, 2011, 50, 9255–9263 CrossRef CAS PubMed.
  4. K. Danyal, D. Mayweather, D. R. Dean, L. C. Seefeldt and B. M. Hoffman, J. Am. Chem. Soc., 2010, 132, 6894–6895 CrossRef CAS PubMed.
  5. Z.-Y. Yang, R. Ledbetter, S. Shaw, N. Pence, M. Tokmina-Lukaszewska, B. Eilers, Q. Guo, N. Pokhrel, V. L. Cash, D. R. Dean, E. Antony, B. Bothner, J. W. Peters and L. C. Seefeldt, Biochemistry, 2016, 55, 3625–3635 CrossRef CAS PubMed.
  6. D. J. Lowe and R. N. F. Thorneley, Biochem. J., 1984, 224, 877–886 CrossRef CAS PubMed.
  7. D. J. Lowe and R. N. F. Thorneley, Biochem. J., 1984, 224, 895–901 CrossRef CAS PubMed.
  8. R. N. F. Thorneley and D. J. Lowe, Biochem. J., 1984, 224, 887–894 CrossRef CAS.
  9. R. N. F. Thorneley and D. J. Lowe, Biochem. J., 1984, 224, 903–909 CrossRef CAS PubMed.
  10. D. F. Harris, Z.-Y. Yang, D. R. Dean, L. C. Seefeldt and B. M. Hoffman, Biochemistry, 2018, 57, 5706–5714 CrossRef CAS PubMed.
  11. P. E. Doan, J. Telser, B. M. Barney, R. Y. Igarashi, D. R. Dean, L. C. Seefeldt and B. M. Hoffman, J. Am. Chem. Soc., 2011, 133, 17329–17340 CrossRef CAS PubMed.
  12. R. Y. Igarashi, M. Laryukhin, P. C. Dos Santos, H.-I. Lee, D. R. Dean, L. C. Seefeldt and B. M. Hoffman, J. Am. Chem. Soc., 2005, 127, 6231–6241 CrossRef CAS PubMed.
  13. D. A. Lukoyanov, N. Khadka, Z.-Y. Yang, D. R. Dean, L. C. Seefeldt and B. M. Hoffman, Inorg. Chem., 2018, 57, 6847–6852 CrossRef CAS PubMed.
  14. J. Christiansen, R. C. Tittsworth, B. J. Hales and S. P. Cramer, J. Am. Chem. Soc., 1995, 117, 10017–10024 CrossRef CAS.
  15. C. Van Stappen, R. Davydov, Z.-Y. Yang, R. Fan, Y. Guo, E. Bill, L. C. Seefeldt, B. M. Hoffman and S. DeBeer, Inorg. Chem., 2019, 58 DOI:10.1021/acs.inorgchem.9b01951.
  16. S. J. Yoo, H. C. Angove, V. Papaefthymiou, B. K. Burgess and E. Münck, J. Am. Chem. Soc., 2000, 122, 4926–4936 CrossRef CAS.
  17. A. D. Scott, V. Pelmenschikov, Y. Guo, L. Yan, H. Wang, S. J. George, C. H. Dapper, W. E. Newton, Y. Yoda, Y. Tanaka and S. P. Cramer, J. Am. Chem. Soc., 2014, 136, 15942–15954 CrossRef CAS PubMed.
  18. L. Cao, O. Caldararu and U. Ryde, J. Chem. Theory Comput., 2018, 14, 6653–6678 CrossRef CAS PubMed.
  19. I. Dance, Inorganics, 2019, 7, 8 CrossRef CAS.
  20. I. Dance, J. Am. Chem. Soc., 2005, 127, 10925–10942 CrossRef CAS PubMed.
  21. I. Dance, Inorg. Chem., 2013, 52, 13068–13077 CrossRef CAS PubMed.
  22. I. Dance, Dalton Trans., 2015, 44, 18167–18186 RSC.
  23. I. Dance, Z. Anorg. Allg. Chem., 2015, 641, 91–99 CrossRef CAS.
  24. B. Benediktsson and R. Bjornsson, Inorg. Chem., 2017, 56, 13417–13429 CrossRef CAS PubMed.
  25. L. Cao, O. Caldararu and U. Ryde, J. Phys. Chem. B, 2017, 121, 8242–8262 CrossRef CAS PubMed.
  26. J. Schlesier, M. Rohde, S. Gerhardt and O. Einsle, J. Am. Chem. Soc., 2016, 138, 239–247 CrossRef CAS PubMed.
  27. O. H. Lowry, N. J. Rosebrough, A. Lewis Farr and R. J. Randall, J. Biol. Chem., 1951, 193, 265–275 CAS.
  28. B. Ravel and M. Newville, J. Synchrotron Radiat., 2005, 12, 537–541 CrossRef CAS PubMed.
  29. T. Spatzal, M. Aksoyoglu, L. Zhang, S. L. A. Andrade, E. Schleicher, S. Weber, D. C. Rees and O. Einsle, Science, 2011, 334, 940 CrossRef CAS PubMed.
  30. S. Metz, J. Kästner, A. A. Sokol, T. W. Keal and P. Sherwood, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2014, 4, 101–110 CAS.
  31. I. T. Todorov, W. Smith, K. Trachenko and M. T. Dove, J. Mater. Chem., 2006, 16, 1911–1918 RSC.
  32. A. D. MacKerell Jr, D. Bashford, M. Bellott, R. L. Dunbrack Jr, J. D. Evanseck, M. J. Field, S. Fischer, J. Gao, H. Guo, S. Ha, D. Joseph-McCarthy, L. Kuchnir, K. Kuczera, F. T. K. Lau, C. Mattos, S. Michnick, T. Ngo, D. T. Nguyen, B. Prodhom, W. E. Reiher III, B. Roux, M. Schlenkrich, J. C. Smith, R. Stote, J. Straub, M. Watanabe, J. Wiórkiewicz-Kuczera, D. Yin and M. Karplus, J. Phys. Chem. B, 1998, 102, 3586–3616 CrossRef PubMed.
  33. A. D. Mackerell Jr, M. Feig and C. L. Brooks III, J. Comput. Chem., 2004, 25, 1400–1415 CrossRef PubMed.
  34. F. Neese, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2018, 8, e1327 Search PubMed.
  35. J. Tao, J. P. Perdew, V. N. Staroverov and G. E. Scuseria, Phys. Rev. Lett., 2003, 91, 146401 CrossRef PubMed.
  36. E. van Lenthe, E. J. Baerends and J. G. Snijders, J. Chem. Phys., 1993, 99, 4597–4610 CrossRef CAS.
  37. C. van Wüllen, J. Chem. Phys., 1998, 109, 392–399 CrossRef.
  38. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305 RSC.
  39. D. A. Pantazis, X.-Y. Chen, C. R. Landis and F. Neese, J. Chem. Theory Comput., 2008, 4, 908–919 CrossRef CAS PubMed.
  40. K. Eichkorn, F. Weigend, O. Treutler and R. Ahlrichs, Theor. Chem. Acc., 1997, 97, 119–124 Search PubMed.
  41. S. Grimme, S. Ehrlich and L. Goerigk, J. Comput. Chem., 2011, 32, 1456–1465 CrossRef CAS PubMed.
  42. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1799 CrossRef CAS PubMed.
  43. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104 CrossRef PubMed.
  44. F. Neese, F. Wennmohs, A. Hansen and U. Becker, Chem. Phys., 2009, 356, 98–109 CrossRef CAS.
  45. S. Kossmann and F. Neese, Chem. Phys. Lett., 2009, 481, 240–243 CrossRef CAS.
  46. R. Izsak and F. Neese, J. Chem. Phys., 2011, 135, 144105 CrossRef PubMed.
  47. T. Petrenko, S. Kossmann and F. Neese, J. Chem. Phys., 2011, 134, 054116 CrossRef PubMed.
  48. R. Bjornsson, F. A. Lima, T. Spatzal, T. Weyhermüller, P. Glatzel, E. Bill, O. Einsle, F. Neese and S. DeBeer, Chem. Sci., 2014, 5, 3096–3103 RSC.
  49. R. Bjornsson, F. Neese and S. DeBeer, Inorg. Chem., 2017, 56, 1470–1477 CrossRef CAS PubMed.
  50. T. Lovell, J. Li, T. Liu, D. A. Case and L. Noodleman, J. Am. Chem. Soc., 2001, 123, 12392–12410 CrossRef CAS PubMed.
  51. O. Einsle, F. A. Tezcan, S. L. A. Andrade, B. Schmid, M. Yoshida, J. B. Howard and D. C. Rees, Science, 2002, 297, 1696–1700 CrossRef CAS PubMed.
  52. K. M. Lancaster, M. Roemelt, P. Ettenhuber, Y. Hu, M. W. Ribbe, F. Neese, U. Bergmann and S. DeBeer, Science, 2011, 334, 974–977 CrossRef CAS PubMed.
  53. S. J. George, R. Y. Igarashi, Y. Xiao, J. A. Hernandez, M. Demuez, D. Zhao, Y. Yoda, P. W. Ludden, L. M. Rubio and S. P. Cramer, J. Am. Chem. Soc., 2008, 130, 5673–5680 CrossRef CAS PubMed.
  54. M. C. Corbett, Y. Hu, F. Naderi, M. W. Ribbe, B. Hedman and K. O. Hodgson, J. Biol. Chem., 2004, 279, 28276–28282 CrossRef CAS PubMed.
  55. S. J. George, B. M. Barney, D. Mitra, R. Y. Igarashi, Y. Guo, D. R. Dean, S. P. Cramer and L. C. Seefeldt, J. Inorg. Biochem., 2012, 112, 85–92 CrossRef CAS PubMed.
  56. L. Rao, X. Xu and C. Adamo, ACS Catal., 2016, 6, 1567–1577 CrossRef CAS.
  57. D. V. Formitchev, C. C. McLauchlan and R. H. Holm, Inorg. Chem., 2002, 41, 958–966 CrossRef PubMed.
  58. P. K. Mascharak, G. C. Papaefthymiou, W. H. Armstrong, S. Foner, R. B. Frankel and R. H. Holm, Inorg. Chem., 1983, 22, 2851–2858 CrossRef CAS.
  59. E. Guan and B. C. Gates, ACS Catal., 2018, 8, 482–487 CrossRef CAS.
  60. D. Yang, S. O. Odoh, T. C. Wang, O. K. Farha, J. T. Hupp, C. J. Cramer, L. Gagliardi and B. C. Gates, J. Am. Chem. Soc., 2015, 137, 7391–7396 CrossRef CAS PubMed.
  61. C.-Y. Fang, S. Zhang, Y. Hu, M. Vasiliu, J. E. Perez-Aguilar, E. T. Conley, D. A. Dixon, C.-Y. Chen and B. C. Gates, ACS Catal., 2019, 9, 3311–3321 CrossRef CAS.
  62. D. Zhao, Z. Chen, W. Yang, S. Liu, X. Zhang, Y. Yu, W.-C. Cheong, L. Zheng, F. Ren, G. Ying, X. Cao, D. Wang, Q. Peng, G. Wang and C. Chen, J. Am. Chem. Soc., 2019, 141, 4086–4093 CrossRef CAS PubMed.
  63. C. Schöttle, E. Guan, A. Okrut, N. A. Grosso-Giordano, A. Palermo, A. Solovyov, B. C. Gates and A. Katz, J. Am. Chem. Soc., 2019, 141, 4010–4015 CrossRef PubMed.
  64. Q.-C. Wang, J.-K. Meng, X.-Y. Yue, Q.-Q. Qiu, Y. Song, X.-J. Wu, Z.-W. Fu, Y.-Y. Xia, Z. Shadike, J. Wu, X.-Q. Yang and Y.-N. Zhou, J. Am. Chem. Soc., 2019, 141, 840–848 CrossRef CAS PubMed.
  65. S. Calvin, XAFS for Everyone, ed. L. Hahn, CRC Press, Taylor & Francis Group, Boca Raton, FL, USA, 2013, vol. 10, p. 281 Search PubMed.
  66. C. Hauser, E. Bill and R. H. Holm, Inorg. Chem., 2002, 41, 1615–1624 CrossRef CAS PubMed.
  67. J. M. Berg, K. O. Hodgson and R. H. Holm, J. Am. Chem. Soc., 1979, 101, 4586–4593 CrossRef CAS.
  68. M. A. Bobrik, E. J. Laskowski, R. W. Johnson, W. O. Gillum, J. M. Berg, K. O. Hodgson and R. H. Holm, Inorg. Chem., 2002, 17, 1402–1410 CrossRef.
  69. G. Moula, T. Matsumoto, M. E. Miehlich, K. Meyer and K. Tatsumi, Angew. Chem., Int. Ed. Engl., 2018, 57, 11594–11597 CrossRef CAS PubMed.
  70. Y. Ohki, K. Tanifuji, N. Yamada, M. Imada, T. Tajima and K. Tatsumi, Proc. Natl. Acad. Sci. U. S. A., 2011, 108, 12635–12640 CrossRef CAS PubMed.
  71. C. R. Sharp, J. S. Duncan and S. C. Lee, Inorg. Chem., 2010, 49, 6697–6705 CrossRef CAS PubMed.
  72. B. B. Wenke, T. Spatzal and D. C. Rees, Angew. Chem., Int. Ed. Engl., 2019, 58, 3894–3897 CrossRef CAS PubMed.
  73. P. Strop, P. M. Takahara, H.-J. Chiu, H. C. Angove, B. K. Burgess and D. C. Rees, Biochemistry, 2001, 40, 651–656 CrossRef CAS PubMed.
  74. M. M. Georgiadis, H. Komiya, P. Chakrabarti, D. Woo, J. J. Kornuc and D. C. Rees, Science, 1992, 257, 1653–1659 CrossRef CAS PubMed.
  75. F. A. Tezcan, J. T. Kaiser, D. Mustafi, M. Y. Walton, J. B. Howard and D. C. Rees, Science, 2005, 309, 1377–1380 CrossRef CAS PubMed.
  76. I. Dance, J. Inorg. Biochem., 2017, 169, 32–43 CrossRef CAS PubMed.
  77. M. C. Durrant, Biochem. J., 2001, 355, 569–576 CrossRef CAS PubMed.
  78. M. J. Dilworth, K. Fisher, C. H. Kim and W. E. Newton, Biochemistry, 1998, 37, 17495–17505 CrossRef CAS PubMed.
  79. L. A. H. van Bergen, M. Alonso, A. Pallo, L. Nilsson, F. De Proft and J. Messens, Sci. Rep., 2016, 6, 30369 CrossRef CAS PubMed.
  80. L. Cao and U. Ryde, Phys. Chem. Chem. Phys., 2019, 21, 2480–2488 RSC.
  81. D. Sippel, M. Rohde, J. Netzer, C. Trncik, J. Gies, K. Grunau, I. Djurdjevic, L. Decamps, S. L. A. Andrade and O. Einsle, Science, 2018, 359, 1484–1489 CrossRef CAS PubMed.
  82. T. Spatzal, K. A. Perez, O. Einsle, J. B. Howard and D. C. Rees, Science, 2014, 345, 1620–1623 CrossRef CAS.
  83. T. Spatzal, K. A. Perez, J. B. Howard and D. C. Rees, eLife, 2015, 4, e11620 CrossRef PubMed.
  84. B. Benediktsson, A. T. Thorhallsson and R. Bjornsson, Chem. Commun., 2018, 54, 7310–7313 RSC.
  85. I. Dance, Dalton Trans., 2016, 45, 14285–14300 RSC.


Electronic supplementary information (ESI) available. See DOI: 10.1039/c9sc02187f

This journal is © The Royal Society of Chemistry 2019