Analysis of early intermediate states of the nitrogenase reaction by Se incorporation and regularization of EPR spectra

doi:10.21203/rs.3.rs-3120611/v1

Download PDF

Article

Analysis of early intermediate states of the nitrogenase reaction by Se incorporation and regularization of EPR spectra

https://doi.org/10.21203/rs.3.rs-3120611/v1

This work is licensed under a CC BY 4.0 License

Journal Publication

published 13 May, 2024

Read the published version in Nature Communications →

Version 1

posted

You are reading this latest preprint version

Due to the complexity of the catalytic FeMo cofactor site in nitrogenases that mediates the reduction of molecular nitrogen to ammonium, mechanistic details of this reaction remain under debate. In this study, selenium- and sulfur-incorporated FeMo cofactors of the catalytic MoFe protein component from Azotobacter vinelandii were prepared under turnover conditions and investigated by using different EPR methods. Complex signal patterns were observed in the continuous wave EPR spectra of selenium-incorporated samples, which were analyzed by Tikhonov regularization, a method that has not yet been applied to high spin systems of transition metal cofactors, and by an already established grid-of-error approach. Both methods yielded similar probability distributions that revealed the presence of at least four other species with different electronic structures in addition to the ground state E₀. Some of these species were preliminary assigned to hydrogenated E₂ states. In addition, advanced pulsed-EPR experiments were utilized to verify the incorporation of sulfur and selenium into the FeMo cofactor, and to assign hyperfine couplings of ³³S and ⁷⁷Se that directly couple to the FeMo cluster. With this analysis, we report selenium incorporation under turnover conditions as a straightforward approach to stabilize and analyze early intermediate states of the FeMo cofactor.

Biological sciences/Biochemistry/Biophysical chemistry

Physical sciences/Chemistry/Physical chemistry/Biophysical chemistry

Nitrogenase

FeMo cofactor

stable isotope labeling

EPR spectroscopy

reaction intermediates

regularization

The conversion of the largely inert N₂ molecule to bioavailable ammonia is essential for life on earth and is a critical step in the biological nitrogen cycle. Biological nitrogen fixation is catalyzed by enzymes of the nitrogenase family that are widespread in bacteria and archaea, but absent in eukaryotes¹. Three isoforms of nitrogenases are distinguished based on the composition of their catalytic cofactor: the Mo-dependent, V-dependent, and Fe-only nitrogenases^2,3. All nitrogenases are two-component proteins consisting of (i) the [4Fe:4S] cluster-containing homodimeric Fe-protein (component Av2 in Azotobacter vinelandii) that serves as reductase and site of ATP hydrolysis, and (ii) the catalytic, α₂β₂-heterotetrameric (or heterohexameric in case of V and Fe) MoFe protein (component Av1 in Azotobacter vinelandii) with two metal cofactors, the [8Fe:7S] P-cluster and the catalytic cofactor. The latter is designated as FeMo cofactor in Mo-dependent nitrogenases and is the most complex bioinorganic metal cluster known to date. The FeMo cofactor consists of seven Fe atoms, nine S atoms, one Mo atom, a central C (carbide) atom, and an organic R-homocitrate moiety (Fig. 1, inset), and is accordingly complex in its electronic and magnetic properties^4–7.

Important traits of the molecular mechanism of nitrogen reduction remain under discussion. It is known that the FeMo cofactor binds the natural substrate N₂ (alternatively also a variety of other small molecules such as CO) during catalysis, and strictly sequentially accepts electrons from the [4Fe:4S] cluster of the Fe protein. This transfer is coupled to the hydrolysis of 2 ATP/e^– by the Fe protein, whereby one electron is first transferred from the reduced P cluster to the FeMo cofactor, and the electron deficit at the P cluster is subsequently replenished by the Fe protein⁸. The reductase component then dissociates from the MoFe protein for reduction and nucleotide exchange before the next 1-electron transfer can take place⁹. Largely due to the complexity of this process, Fe protein is the only known reductant to sustain productive N₂ reduction by MoFe protein, although recent electrochemical approaches have been reported to achieve similar results¹⁰. The reduction of N₂ follows a minimal stoichiometry of

N₂ + 8 e^– + 8 H⁺ + 16 ATP ࣧ→ 2 NH₃ + H₂ + 16 [ADP + P_i],

including the obligatory release of H₂ with a limiting stoichiometry of 1 H₂/N₂. The kinetics of the reaction are comprehensively outlined in a scheme proposed by Lowe and Thorneley (LT)¹¹, in which the system cycles through eight distinct states, E₀ to E₇, each representing the addition of a single electron (Fig. 1). Under reductive conditions the FeMo cofactor is commonly isolated in the resting state E₀¹¹, and then successively receives electrons (and protons for charge balance) through states E₁ to E₇. Importantly, the binding and activation of N₂ requires the enzyme to reach state E₃ or E₄, which is complicated by the risk of an unproductive loss of 2 electrons as additional H₂^2,12. This finding indicated that an essential aspect of electron accumulation on the FeMo cofactor is the formation of surface hydrides that can be lost as H₂ by accidental protonation^3,13. Stabilization of these surface-associated hydride adducts may be achieved by a bridging binding mode¹⁴; this type of electron storage is crucial for the cluster to accumulate four electrons at isopotential (i.e., from the Fe protein) and allows for a mechanistic twist upon reaching the E₃ or E₄ state. Triggered by the presence or binding of the substrate N₂, the two adjacent hydrides present in the E₄ state can reductively eliminate H₂, leaving the enzyme in a 2-electron-reduced state that cannot be achieved by electron transfer from the Fe protein alone and that is sufficiently reactive to break the N₂ triple bond¹⁵. From states E₅-E₇, the reaction then proceeds to the release of the product NH₃, but different mechanistic routes remain under debate^2,16,17.

--- Please insert Fig. 1 here ---

The E₀ state of the FeMo cofactor has a total spin of S = 3/2^18,19 and the oxidation state of the FeMo cofactor changes by 1 with each reaction step; so that the total spin of the cofactor is half-integer for any even state and integer for any odd state. The odd states are thus either diamagnetic (S = 0) or have "non-Kramers" spin states² with high zero-field splitting and hence the absence of EPR transitions at common EPR frequencies^20,21. EPR spectroscopy provides access to the characterization of the ground state as well as the S\(\ne\) 0 reaction intermediates, and supports the drawing of mechanistic and – within limits – also structural conclusions. In particular, freeze-quenched samples with different substrates, some of them stable-isotope-labeled, have been studied^22–24. Several of these studies showed complex continuous-wave (cw)-EPR spectra with well-resolved anisotropy of the g-tensor, indicating that the substrate directly couples with at least one Fe atom of the FeMo cofactor^24,25. However, an unambiguous assignment of the binding position was not possible.

The substrate binding site of the CO-inhibited FeMo cofactor in its resting state was identified by crystallography^26,27. CO displaces the S at position S2B and a CO bond in an end-on µ₂-bridging mode to Fe2 and Fe6 is formed at this position^26,27. In a subsequent study, KSeCN was found to be both a substrate and an inhibitor of nitrogenase activity, and crystal structures from freeze-quenched nitrogenase samples generated during turnover with KSeCN revealed that S2B was replaced by Se²⁸. When KSeCN was removed from the reaction mixture and the reaction was allowed to proceed, further Se exchange first occurred at positions 3A and 5A, the other two µ₂-bridging S that form the equatorial ‘belt’ of the cofactor (Fig. 1, inset). Only after several thousand more reaction cycles, the incorporated Se was again replaced by S. Starting from the exclusive Se2B labeling, the approximately equal labeling distribution of the other two positions (3A and 5A) was reached after about 1000 turnover cycles²⁸. Comparable S-to-S exchange experiments within the sulfur belt were carried out with the VFe cofactor of V-dependent nitrogenase²⁹. A subsequent study examining Se incorporation into the FeMo cofactor of a Mo-dependent nitrogenase at high and low KSeCN concentrations established that both conditions lead to a similar Se distribution within the cofactor³⁰. Furthermore, it could be demonstrated that Se labeling is also possible at positions 3A and 5A by gassing the 2B-Se-labeled protein with CO during catalysis. In this process, the Se2B is exchanged by CO, while the two S atoms at the 3A/5A positions are replaced by Se. The use of such a Se-labeled FeMo cofactor allowed its electronic structure to be analyzed by various methods like X-ray spectroscopy³⁰.

Based on these studies, the goal of this work was to determine whether and to what extent Se is incorporated into the FeMo cofactor and what geometric or electronic changes result from this manipulation. We use high-resolution EPR spectroscopy for this purpose, as structure-determination methods can identify the labeling positions of individual isotopes within the FeMo cofactor, but the various electronic structures or redox states of the cluster are difficult to be distinguished other than with complex approaches like spatially resolved anomalous dispersion refinement³¹. Tikhonov regularization, commonly applied to analyze complex magnetic resonance datasets, e. g., from PELDOR/DEER spectroscopy^32–36, was employed for the first time on cw-EPR spectra of the high-spin FeMo cofactor to assign individual species formed by Se incorporation. The resulting probability distributions revealed several species with different electronic structures in each sample, making an assignment to specific intermediates and/or redox states possible. The quality of our analyses was compared to those obtained from a grid-of-error approach³⁷.

Together, these studies establish that Se incorporation into the FeMo cofactor provides access to other states in the kinetic LT scheme that will help to better understand the molecular mechanism of the FeMo cofactor in the nitrogenase reaction.

Sample preparation. The MoFe-protein and Fe-protein from Azotobacter vinelandii (designated as Av1 and Av2, respectively) were isolated under anoxic conditions as described previously²⁸.

Enzyme assays. Turnover assays for Av1 and Av2 were prepared in a buffer containing 50 mM Tris-HCl (pH 7.5), 200 mM NaCl, 5 mM Na₂S₂O₄ and supplemented with 20 mM creatine phosphate, 5 mM ATP, 5 mM MgCl₂, 25 units/mL phosphocreatine kinase and 25 mM Na₂S₂O₄ (in 50 mM Tris-HCl, pH 7.5 and 200 mM NaCl)²⁸. All samples except for the Av1-Se-C₂H₂ sample were kept under an argon/H₂ atmosphere and the indicated amounts of KSCN, K³³SCN, KSeCN, or K⁷⁷SeCN were added to the reaction (see Table 1). C₂H₂ was used as substrate in the Av1-Se-C₂H₂ sample. Afterward, the Av2 protein and remaining SCN^– or SeCN^– were removed by three rounds of sample concentration and dilution with a 100-kDa molecular weight cut-off ultrafiltration device (Vivaspin, Sartorius). An additional desalting step (Sephadex G25, GE Healthcare) was applied with samples Av1-WT, Av1-³³S, Av1-Se2B-1, Av1-Se-C₂H₂, Av1-Se-low and Av1-⁷⁷Se2B. Sample concentrations were determined by absorbance at 410 nm²⁸; relative EPR signal intensities were determined by double-integration of the respective X-band cw-EPR spectra.

Cw-EPR experiments. X-band cw-EPR experiments were performed using Bruker E500 or E580 spectrometers in combination with Bruker resonators (4122SHQE or 4119HS-W1) combined with an Oxford ESR900 helium gas flow cryostat. Power-sweep experiments were done at 5 K, a microwave frequency of 9.39 GHz, a modulation amplitude of 0.6 mT, and a conversion time of 165.25 ms. For testing the relaxation behavior of the individual samples, cw-EPR spectra at different microwave powers (from 0.025 to 39.4 mW at the E500, or from 0.377 to 37.7 mW at the E580) were recorded. Temperature-dependent experiments were recorded at 6, 9, or 12 K using a microwave power of 0.095 mW, a modulation amplitude of 0.6 mT, and a conversion time of 165.25 ms.

Light induced cw-EPR experiments. Similar to the protocol described in¹⁴, two samples, Av1-Se2B-1 and Av1-Se-low, were illuminated inside the cooled cavity (Bruker 4119HS-W1) in combination with the cryostat (Oxford ESR900) for about 10 min using a blue-light LED (100 mW, Schott KL 2500). The cw-EPR experiments were performed at 6 K at 9.38 GHz by using microwave power of 3.77 mW, a conversion time of 160 ms and a modulation amplitude of 0.6 mT. The cryogen annealing was done by keeping the samples in a cryogen-solution (isopropanol-liquid nitrogen) at about 150 K for some hours. Additionally, the samples were stored for 16 h in liquid nitrogen.

Pulse EPR experiments. Pulse Q-Band EPR experiments were performed using a Bruker E580 spectrometer in combination with a Bruker EN 5107D2-flexline resonator immersed in an Oxford CF935 helium gas-flow cryostat. All experiments were carried out at a microwave frequency of 33.8 GHz at 4.5 K. Unless noted otherwise, a video gain setting of 200 MHz was used.

Longitudinal transient nutation experiments. Experimental conditions: pulse length π/2 = 10 ns, nutation step width 4 ns, τ = 110 ns, T = 600 ns, and a shot repetition time of 51 µs. A 4-step phase cycle was used. The spectra were measured in steps of 10 mT. As the nutation frequency depends on the local microwave magnetic field strength B₁, all frequency axes were normalized to the nutation frequency measured with a coal reference sample (Bruker). This standardization makes the frequency axis essentially independent of spectrometer-specific settings such as microwave power or the resonator quality (Q-factor). The nutation signals were processed as follows: After subtraction of a polynomial baseline, a Hamming window function and a zero filling with a fill factor of 4 were applied. Finally, an FFT was performed.

Inversion recovery experiments. Experimental conditions: pulse lengths π/2 = 12 ns, τ = 100 ns, T_start = 400 ns, T-steps = 80 ns, and a shot repetition time of 100 µs. The video gain was set to 20 MHz. The spectra were measured in steps of 3 mT. From each spectrum the resonator background was subtracted. Exponential fit functions were used to determine \({T}_{1}^{\text{e}\text{f}\text{f}}\).

2-pulse ESEEM versus B ₀ experiments. Experimental conditions: pulse length π/2 = 12 ns, τ_start = 100 ns, τ-steps = 4 ns with 40 steps and a shot repetition time of 20 µs. The spectra were measured in steps of 0.3253 mT. The resonator background was subtracted from each spectrum. Pseudo modulation was performed using a modulation amplitude of 1.0 mT and a binominal smoothing with 4 smoothing points.

For determining \({T}_{\text{M}}^{\text{e}\text{f}\text{f}}\), modified experimental conditions were used: pulse length π/2 = 12 ns, τ_start = 100 ns, τ-steps = 4 ns with 500 steps and a shot repetition time of 50 µs. A 16-step phase cycle was used. The spectra were measured in steps of 3.0 mT. Exponential fit functions were used for analysis.

3-pulse ESEEM experiments. Experimental conditions: pulse lengths π/2 = 10 ns, τ = 90 ns, T_start = 100 ns, T-steps = 8 ns with 750 steps, shot repetition time 70 µs. The spectra were measured in steps of 10 mT. A 4-step phase cycling was used. Spectra have been processed as follows: The phase of the time domains have been optimized, a mono- or bi-exponential background function has been subtracted, a Hamming window function has been applied, a zero-filling factor of 4 has been used, and finally, a cross-term-averaged FFT was applied.

Data Analysis. Spectral simulations of cw-EPR spectra were carried out using the Matlab (The MathWorks, Natick, MA) package EasySpin ³⁸ with its “pepper” simulation routine; spectral analysis was done using self-written Matlab scripts. The regularization and the grid-of-errors method were implemented as Matlab scripts (for details see Supporting Information, part A). The regularization results were analyzed using a multi-Gaussian approach described in the Supporting Information, section B11.

3P-ESEEM simulations were carried out using the EasySpin algorithm “saffron”³⁹. Pseudo-nuclear and effective hyperfine couplings were included by calculating with a total electron spin quantum number of S = 3/2 and zero-field coupling D = 180 GHz. ESEEM signals of the two nitrogen atoms were simulated using literature parameters: A(N1) = [1.02 0.98 1.14] MHz, Q(N1) = 2.17 MHz/η(N1) = 0.6, and A(N2) = [0.5 0.4 0.4] MHz, Q(N2) = 3.5 MHz/η(N2) = 0.35; Euler angles of 60°, 20°, 0° between the g and quadrupolar tensor for the second nucleus axis were used⁴⁰.

Spectral simulations of ESEEM signals of sample Av1-⁷⁷Se2B were performed as follows: Using the determined ¹⁴N hyperfine couplings and assuming one ⁷⁷Se nucleus, the analysis of the spectral pattern in the Av1-⁷⁷Se2B sample was done by manual optimization. For a one-to-one simulation, different rhombicity (\(\lambda\)) values that affect both the effective hyperfine couplings and the pseudonuclear g-factors were taken into account. These effects scale with the magnitude of the hyperfine couplings and thus, alter the effective nuclear Larmor frequency and the effective hyperfine couplings. For this reason, simulations were performed for each \(\lambda\) value individually. The simulations were done between 0 ≤ \(\lambda\) ≤ 1/3 in 167 steps. For each \(\lambda\) value the 3P-ESEEM spectra S(\(\lambda\)) were calculated and weighted by the probability-distribution P(\(\lambda\)) obtained from regularization. The total spectrum was obtained by: \(S={\sum }_{\lambda }S\left(\lambda \right)P\left(\lambda \right)\). Only the lower Kramers doublet was considered.

Runtime estimation of the regularization and grid-of-errors methods. Using a standard desktop PC (Intel Core i5-4590 CPU @ 3.3 GHz) with Matlab 2019a and EasySpin 5.2.25, the calculation of the kernels (667 \(\lambda\)-steps with \(0\le \lambda \le 1/3\), 18 intrinsic lineshape-steps with 0.5 mT steps, and an angle step width of 0.5° for calculation of a powder spectrum) required about 12 hours. The regularization itself, using 27 α-values and 18 \(lwpp\) points required a compute time of approximately 150 seconds. It is therefore time-saving to calculate the kernel once per series of spectra. On the other hand, the calculation of the grid (223 \(\lambda\)-steps with \(0\le \lambda \le 1/3\), 249 \(lwpp\)-steps \(0 \text{m}\text{T}\le lwpp\le 25 \text{m}\text{T}\), and an angle step width of 0.5° for calculation of a powder spectrum) required about 56 hours. The grid-of-errors optimization itself required only about 300 seconds. The comparison of the compute times clearly shows that the regularization method requires less computing time and should therefore be preferred over the grid-of-errors method if the prerequisites for regularization are fulfilled (see below). The reduction of computation time is mainly due to the lower required number of steps in the second parameter dimension (here: \(lwpp\)). By choosing identical number of steps for \(\lambda\) and \(lwpp\), the compute times for both methods are very similar.

Regularization of cw-EPR spectra. To spectroscopically follow the changes of the FeMo cofactor after incorporation of Se, different nitrogenase Av1 samples were produced under various turnover conditions in the absence of N₂ (see Table 1); these samples exhibited different labeling positions (position 2B and/or positions 2B, 3A, and 5A) or labeling yields. For this purpose, different KSeCN and KSCN concentrations (samples Av1-Se2B-1, Av1-Se-low, Av1-Se2B- lowflux), different Av1/Av2 ratios (samples Av1-Se2B-lowflux and Av1-S) and different reaction cycles (samples Av1-Se-C₂H₂ and Av1-S-remigration) were applied. Two S-incorporated samples, one with ³³S (Av1-³³S) and one with natural abundance ³²S (Av1-S) were prepared under turnover conditions and analyzed in comparison. All samples were frozen after the defined number of reaction cycles, but not under freeze-quench conditions. Therefore, no short-lived intermediate states are expected to be trapped. Figure 2 depicts the cw-EPR spectra of all Av1 samples under investigation covering a magnetic field range of 50–283 mT.

--- Please insert Fig. 2 here ---

The Av1-WT sample in its resting state exhibits the well-known S = 3/2 spin state EPR spectrum of the lower Kramer’s doublet (panel A). The two EPR spectra of the S-incorporated samples, Av1-S and Avl-³³S (panels B and C) are virtually identical compared to the unmodified protein; therefore, incorporation of S and in particular ³³S (with a nuclear spin of I = 3/2) into the FeMo cofactor is not detectable by cw-EPR spectroscopy. All Se-exchanged samples, however, exhibit a complex signal shape with at least five peaks spanning the 120–260 mT magnetic field range. Unexpectedly, the “Se-patterns” of samples Av1-Se2B-1, Av1-Se2B-lowflux, Av1-Se-C₂H₂, Av1-Se-low, and even of Av1-⁷⁷Se2B (panels D-H) are similar, only differences in individual peaks intensities can be observed. It is important to note that ⁷⁷Se has a nuclear spin of I = ½, which is different to the I = 0 of the naturally most abundant isotopes ⁷⁸Se and ⁸⁰Se. As those samples show very similar spectral patterns, hyperfine couplings of ⁷⁷Se and the FeMo cofactor can be excluded as the origin of the Se-pattern. The cw-EPR spectrum of the Av1-S-remigration sample (panel I) again exhibits the Se-pattern, but with decreased intensity. Qualitatively, the observed signal pattern can be described as a mixture of signals from unlabeled and Se-incorporated samples.

For a more quantitative evaluation of S-, Se- and unlabeled samples, the intensity differences of the respective cw-EPR spectra were compared using spin counting via double integration. Samples Av1-WT, Av1Se2B-1, Av1-⁷⁷Se2B, Av1-Se-low, Av1-Se-C₂H_2, and Av1-³³S were compared, as all were prepared from the same enzyme batch and under identical electron flux. The analysis shows that the signal intensity of sample Av1-³³S is comparable to the intensity of the Av1-WT sample, but all Se-incorporated samples have only ≈ 60% of the resting-state intensity (Supporting Figure B2). Consequently, Se incorporation leads to ≈ 40% EPR-inactive (S = 0) and/or non-Kramers states (S = 1, 2, 3, ...).

It is essential to know the origin of the complex Se-pattern to perform correct spectral simulations of the experimental data. Hyperfine couplings have already been ruled out as the source, geometric distortions due to Se incorporation are also unlikely as the only explanation, as there is no evidence for such in the crystal structures²⁸, assuming that the Se incorporation in crystals is representative of that in solutions. Moreover, the EPR signal pattern of sample Av1-Se-C₂H₂, in which Se should be incorporated over the entire sulfur belt, is almost identical to those of the other Se-incorporated samples labeled mainly at the 2B position (see also below). Therefore, different states of the FeMo cofactor that manifest in different zero-field splitting parameters are the most plausible assumption. In this case, the cw-EPR spectra of all samples are dominated only by the rhombicity parameter (\(\lambda )\) of the zero-field splitting as the effective g-factors \({g}_{\left\{x,y,z\right\}}^{1/2}\) of the lower Kramer doublet of an S = 3/2 system are functions of \(\lambda =|E/D|\) (see Supporting Information Part A, Theory).

Exact \(|E/D|\) values are thus desired for a precise simulation of pulsed EPR data as the zero-field Hamiltonian \({H}_{\text{Z}\text{F}\text{S}}\) depends on \(\left|D\right|\) and \(\lambda =|E/D|\). \(\left|D\right|\) can be estimated experimentally by temperature-dependent measurements of the intensity ratios of the lower and upper Kramers doublet at g ≈ 6⁴¹. These measurements were conducted on samples Av1-WT and Av1-Se2B-1 at 6 K and 15 K (Supporting Figure B3), and small differences were observed: The signal of the latter sample is slightly shifted to ≈ 115 mT and shows a more complex signal pattern compared to the single signal at 111 mT in the Av1-WT sample. However, quantitative extraction of signal intensities was not possible due to the substantial overlap of the signals from the lower and upper Kramer doublet (Supporting Figure B3). Nevertheless, the analysis demonstrates that \(\left|D\right|\) is of the same magnitude in the Se-incorporated samples and hence, using the WT value of \(D=180 \text{M}\text{H}\text{z}\) is a valid approximation. Please note that the effective \(g\)-values are independent of \(D\), if the energy of the Zeeman interaction is small compared to zero-field energy.

Inhomogeneous broadening of the magnetic parameters of protein-bound (metal) cofactors is usually approximated by a random distribution of the EPR parameters, in particular the \(\varvec{D}\) tensor and the \(\varvec{g}\) matrix, using Gaussian distributions, so-called strain models^42–45. These distribution models are valid as long as the width of the distribution is small compared to its magnitude. However, the experimental spectra of the high-spin Se-FeMo cofactor exhibit a large splitting compared to their size (Fig. 2), so that such simple strain models cannot correctly reproduce these data sets, and thus other approaches are required.

Having only the parameter \(\lambda\) that dominates the cw-EPR spectrum, a regularization method was applied to disentangle the complex signal pattern in the Se-incorporated samples (see Supporting Information Part A for theoretical details). Briefly, ill-posed problems can be solved by Tikhonov regularization, and although this method is commonly used in the analysis of DEER datasets^33,34, its application to cw-EPR spectra of high-spin transition metal clusters is yet not established. First, the potential and robustness of the regularization method was thoroughly tested using three calculated model datasets (Supporting Table 1). After the optimal regularization parameter \({\alpha }_{\text{O}\text{p}\text{t}}\) was determined by different methods, the distribution function was obtained. From this, the respective cw-EPR spectrum was calculated (Supporting Figures A3–12). The regularization reproduced the calculated model spectra very well (Supporting Figure A9–12), and therefore, the method was used to analyze all experimental Av1 cw-EPR spectra. As the regularization allows only one free parameter (\(\lambda\)), an intrinsic linewidth (lwpp) analysis of all samples was first performed, and optimal intrinsic Lorentzian peak-to-peak line shapes of 2.5–3 mT, 3.0–3.5 mT, and 3.5–4.0 mT were obtained for spectra recorded at 5–6 K, 9 K and 12 K, respectively (Supporting Information Part A and Supporting Figures B4–8). The distribution functions obtained from regularization are shown in Fig. 3 and the individual \(\lambda\) values of all species are summarized in Table 2. A multi-Gaussian fit was applied to quantify the individual distributions (Supporting Figure B11 and Table B1).

It is observed that samples Av1-WT, Av1-S, and Av1-³³S (panels A–C) contain only one spin species with an average value of \({\lambda }_{2}\) = 0.054. In contrast, samples Av1-Se2B-1, Av1-Se2B-lowflux, Av1-Se-C₂H₂, Av1-Se-low and Av1-⁷⁷Se2B (panels D-H) contain five species with average \(\lambda\) values of \({\lambda }_{1}\) = 0.033, \({\lambda }_{2}\) = 0.057, \({\lambda }_{3}\) = 0.082, \({\lambda }_{4}\) = 0.116 and \({\lambda }_{5}\)≈ 0.19. The second value, \({\lambda }_{2}\), matches that of the Av1-WT and accordingly was assigned to the electronic resting state of the FeMo cofactor. Even though the other four "Se-species" are present in all Se-incorporated samples, noticeable population differences between samples can be detected. In Av1-Se2B-1 and Av1-⁷⁷Se2B, all four Se-species are populated, with \({\lambda }_{4}\) being the largest fraction (~ 36%). In Av1-Se-low, on the other hand, the fraction of species \({\lambda }_{2}\) is below 10%, the Se-species are more highly populated, in particular \({\lambda }_{4}\). It is worth noting that the \(\lambda\) populations of samples Av1-Se2B-1 and Av1-Se2B-lowflux differ; In contrast to Av1-Se2B-1, sample Av1-Se2B-lowflux shows predominantly \({\lambda }_{2}\) and only small amounts of any of the Se-species. This can be rationalized by a lower electron flux in sample Av1-Se2B-lowflux due to the lower Av2/Av1 ratio, which in turn might result in a decreased formation rate of Se-species per time. The largest \({\lambda }_{5}\) value of ≈ 0.19 has a very broad \(\lambda\) distribution and in most cases only a low (< 10%) population. Sample Av1-S-remigration (panel I), in which the Se is expected to be re-replaced by S, shows a different distribution than any of the other Se-incorporated samples: Species \({\lambda }_{1}\) and \({\lambda }_{4}\) are depopulated, and in addition to the resting state, only the \({\lambda }_{3}\) state is populated.

To evaluate the relaxation behavior of the individual spin species, cw-EPR spectra were recorded at different microwave powers of 0.377 mW, 3.77 mW, and 37.7 mW for analysis by regularization (Fig. 3, red and blue lines, additional microwave powers are shown as Supporting Figure B9). The relaxation behavior of all Se-species is similar, but different from that of the resting-state FeMo cofactor (\({\lambda }_{2}\)). Temperature-dependent measurements at 6, 9, and 12 K produced similar results (Supporting Figure B5–8).

--- Please insert Fig. 3 here ---
--- Please insert Table 2 here ---

From the normalized population distributions (Fig. 3), cw-EPR spectra were calculated (red lines in Fig. 2). The agreement between experiment and regularization is remarkably good in all samples and demonstrates the potential of the regularization method. Slight differences, e.g., in the signals at 145 mT and 200 mT (panels D-I), are only intensity differences and are most likely caused by small baseline artifacts.

Analysis of cw-EPR spectra using the grid-of-error method. The question remains whether the cw-EPR spectra are dominated only by the \(\lambda\) parameter or whether the intrinsic line shape lwpp is a second important parameter that differs between samples and/or between individual spin species. Therefore, the established grid-of-error approach³⁷ was used as a second method to re-evaluate all Av1 cw-EPR spectra. The results are depicted in Fig. 4 and demonstrate that this method yields similar distribution functions compared to the regularization method. It is noteworthy that the P(\(\lambda\)) functions are significantly narrower than those obtained by regularization. This is not surprising, as the width of the distribution is partially compensated by a distribution of the intrinsic spectral linewidths. Again, samples Av1-WT Av1-S and Av1-³³S (panel A–C) contain only one species with a \(\lambda\) = 0.054 value, and samples Av1-Se2B-1, Av1-Se2B-lowflux, Av1-Se-C₂H₂, Av1-Se-low and Av1-⁷⁷Se2B (panels D–H) contain four Se-species with \(\lambda\) values of \({\lambda }_{1}\) = 0.035, \({\lambda }_{2}\) = 0.058, \({\lambda }_{3}\) = 0.085, \({\lambda }_{4}\) = 0.12. A fifth species with a \(\lambda\) value of around ≈ 0.19 can be detected in samples Av1-Se2B-1, Av1-Se-C₂H₂, Av1-Se-low, and Av1-⁷⁷Se2B. Sample Av1-S-remigration (panel I) shows only three species with \(\lambda\) values of 0.058, 0.085, and 0.12. These \(\lambda\) values are very similar to those obtained by regularization.

Qualitatively, both methods yield similar population trends for all Se-incorporated samples. However, the individual populations differ depending on the method of analysis, and as we believe that the regularization provides more reliable populations, only for this method, a quantitative evaluation was carried out (Table 2). One major advantage of the grid-of-error method is that two (or even more) parameters can be optimized simultaneously so that linewidths are obtained for all species analyzed. A 2-dimensional representation (\(\lambda\) and lwpp) shows that the non-Se-incorporated cofactors exhibit a lwpp between 1 mT and 3 mT (Supporting Figure B10), consistent with the result of 2.5 mT from regularization. The analysis of the Se-incorporated samples confirms that the lwpp of \({\lambda }_{1}\), \({\lambda }_{2}\) and \({\lambda }_{3}\) are between 1–3 mT, and only the lwpp of \({\lambda }_{4}\) is significantly larger than 5 mT. This result is unexpected, as the analyses of the relaxation times led to similar values for all Se-incorporated samples (see below). One explanation could be that the bandwidth of the individual \({\lambda }_{4}\) values is significantly broader than \({\lambda }_{1-3}\), mainly because the grid-of-error method tends to overrate the parameter lwpp (see also section “Regularization versus grid-of-error approach”).

--- Please insert Fig. 4 here ---

Light excited experiments.

Hoffman and coworkers¹⁴ have used intra-EPR cavity photolysis at 450 nm to characterize hydride containing states of the FeMo cofactor; by irradiating nitrogenase samples with blue light and subsequent annealing at 150 K, a conversion of two E₂(2H) isomers (denoted as 1b and 1c) could be demonstrated. Following these studies, samples Av1-Se2B-1 and Av1-Se-low were used to perform such experiments. The respective cw-EPR spectra (Supporting Figure B12) were analyzed by regularization and are shown in Fig. 5. It is evident that both samples respond to light irradiation and subsequent cryo-annealing, i.e. the probability distributions of the species change, but the changes are more pronounced in sample Av1-Se-low. This may be due to the fact that this sample contains a higher Se concentration.

In contrast to the results presented in reference¹⁴, no species appear upon light illumination, but rather only reduction of signal intensities can be detected (blue arrows in Fig. 5). A one-to-one correspondence to the published results cannot be expected, however, as the FeMo cofactor used in reference¹⁴ and the Se-FeMo cofactors and accompanying intermediates studied in our experiments do have slightly different properties such as binding strengths and absorption coefficients. The regularization clearly shows that the population probabilities of the individual species are different: while \({\lambda }_{2}\) and \({\lambda }_{3}\) do not change, the population probabilities of \({\lambda }_{1}\) and \({\lambda }_{4}\) decrease significantly, and similarly. As the ground state \({\lambda }_{2}\) is not supposed to change, we can identify two distinct responses: The population probabilities of \({\lambda }_{1}\) and \({\lambda }_{4}\) change with light, those of \({\lambda }_{3}\) do not.

--- Please insert Fig. 5 here ---

Pulse EPR experiments. Prior analyses of hyperfine couplings, transient nutation, inversion recovery, and 2-pulse ESEEM experiments were conducted at Q-band microwave frequencies to determine the relaxation times and spin states of all samples. The transient nutation experiments revealed that unlabeled and Se-incorporated samples contain the same nutation frequencies, and only the intensities and linewidths of individual signals differ to a small extent (Supporting Figure C1). Therefore, all “Se-species” must possess the same total spin as the FeMo cofactor in its resting state (S = 3/2). Analysis of 2-pulse ESEEM and inversion recovery spectra yielded the relaxation times \({T}_{\text{M}}^{\text{e}\text{f}\text{f}}\) and \({T}_{1}^{\text{e}\text{f}\text{f}}\), which are in the range of 200–400 ns and 1–3 µs, respectively (Supporting Figures C2 and C3). The relaxation times of all samples are similar, and are too short to conduct certain pulse experiments like ENDOR spectroscopy under our experimental conditions.

Representative Q-Band \(\tau\)-averaged 2-pulse ESEEM experiments of samples Av1-WT and Av1-Se2B-1 are depicted as upper panels of Fig. 6. Additionally, the pseudo-modulated spectra are shown for a direct comparison with the cw-EPR spectra shown in Fig. 2A/D. Both spectra are quite similar to the ones obtained from X-band microwave frequencies: the Av1-WT sample shows the typical spectrum of the FeMo cofactor in its resting state (Fig. 6, left), and the Av1-Se2B-1 sample shows the already described complex Se-pattern (Fig. 6, right). However, the signal-to-noise ratio (S/N) of the pulse EPR spectrum is significantly lower, which is mainly due to the lock-in detection of the cw-EPR spectra, and the intensities of the individual signals differ slightly due to the incomplete compensation of different ESEEM modulation depths at different magnetic field positions by \(\tau\)-averaging.

3P-ESEEM spectra (Fig. 6, lower panels) of Av1-WT (black traces), Av1-Se2B-1 (red traces), and Av1-S (dark blue traces) are depicted at four different magnetic-field positions (580, 660, 740, and 880 mT), these spectra show nearly identical hyperfine couplings close to the proton Larmor frequency and in the range between 0–5 MHz; the latter signals have been assigned to two nitrogen atoms of the surrounding amino acids^40,46. Using literature values^40,46, the ESEEM signals of the three samples can be simulated with good agreement. This result confirms that the direct protein environment of the FeMo cofactor remains structurally intact after turnover with KSeCN, and that no other ligand such as SeCN^– or CN^– is attached to the cluster. In addition, it is reconfirmed that the overall spin of the cluster remains the same, otherwise, additional nitrogen hyperfine couplings would be expected.

On the other hand, samples Av1-⁷⁷Se2B (orange traces) and Av1-³³S (light blue traces) show additional resonances (shaded orange and light blue areas in Fig. 6), which originate from hyperfine couplings of the respective EPR-active nuclei (³³S and ⁷⁷Se) and the FeMo cofactor. Differences in the frequencies and signal patterns are due to different Larmor frequencies of the two nuclei and additional quadrupole couplings in the case of ³³S. Sample Av1-Se2B-1 does not show any Se hyperfine couplings as the natural abundance of ⁷⁷Se is below 8%. Spectral simulations of these additional hyperfine couplings are required for a quantitative analysis. However, such simulations are complex because at almost all magnetic positions the EPR spectra of the Se-species overlap, and therefore the observed ³³S and ⁷⁷Se hyperfine couplings are the weighted sum of each species’ contribution.

Additional difficulties arise when simulating the ³³S hyperfine couplings in sample Av1-³³S, as the quadrupole coupling of the ³³S nucleus overlaps strongly with the resonances of the two ¹⁴N nuclei. Moreover, depending on the magnetic-field position, different ESEEM resonances are suppressed due to cross-suppression effects, and the 3P-ESEEM spectrum of two ¹⁴N nuclei and one ³³S nucleus shows a large number of peaks due to the product rule. Therefore, no unequivocal spectral simulation could be achieved. Qualitatively, the few signals in the 580 mT and 660 mT spectra indicate that a single ³³S nucleus with hyperfine and quadrupole couplings of a few MHz can generate such a pattern.

Using published ¹⁴N hyperfine couplings and assuming one ⁷⁷Se nucleus, the analysis of the spectral pattern in the Av1-⁷⁷Se2B sample was done by manual optimization (see Methods section for details) and yielded principal ⁷⁷Se hyperfine coupling values of A_x = 3 MHz, A_y = 10.5 MHz and A_z = 0 MHz (a_iso (⁷⁷Se) ~ 4 MHz) (grey shaded dotted traces in Fig. 5). Of these values, only A_y can be trusted, as B₀ = 560 mT corresponds to the effective g_y principal value of the \({\lambda }_{2}\) species. Variations of A_x and A_z, especially at higher magnetic fields, do not affect the quality of the simulations, so both values are undefined.

--- Please insert Fig. 6 here ---

Regularization versus grid-of-error approach. For the analysis of the complex cw-EPR spectra, two model-free methods, the grid-of-errors method³⁷ and the regularization method, were chosen to identify and analyze the individual spin species. The former method has been successfully applied to a high-spin Fe-EDTA complex^37,47. An accurate \(|E/D|\) value is necessary for both methods, as only then the computed rhombicity values can be converted into a correct effective g-matrix (see Supporting Information Part Theory). In the Fe-EDTA system, \(\left|D\right|\) is not significantly larger than the electron-Zeeman splitting in X-band, therefore measurements at several magnetic field strengths and simultaneous evaluation of all spectra with the grid-of-errors approach lead to accurate \(\left|D\right|\) values⁴⁷. In the FeMo cofactor, \(\left|D\right|\) (\(\approx\) 180 GHz¹⁹) is much larger than the electron Zeeman splitting in X-band (\(\approx 10 \text{G}\text{H}\text{z})\) and therefore, \(\left|D\right|\) can only be precisely determined at frequencies above \(\left|2D\right|\), or by performing a frequency sweep experiment at different magnetic fields⁴⁸. Such experiments are quite difficult to perform in terms of sample size and experimental conditions; however, the qualitative analysis performed in this study showed that \(\left|D\right|\) can be safely assumed unchanged in all samples (Supporting Figure B3).

As regularization has never been applied to statistical distributions of the zero-field parameter in high-spin systems, both analytical methods were first tested and compared using three model systems. A fixed intrinsic lineshape lwpp of 1 mT was used, thus the only variable parameter in these simulations was the rhombicity \(\lambda\). Comparison of the calculated and simulated spectra showed that the grid-of-errors approach, in particular in the case of low S/N, gave inferior results in comparison to the regularization method, which consistently performed exceedingly well (Supporting Information Part A).

For analysis of the experimental FeMo cofactor spectra, a lwpp of 2.5 mT was determined for all samples at 5 K from a lwpp analysis (determination of the minimum in a \(\rho \left(lwpp\right)\) versus \(lwpp\) plot) and was used in the regularization method (Supporting Figure B4). The lwpp was used as a second independent parameter in the grid-of-error approach; this may be advantageous if the lwpp differs from species to species. In the Se-incorporated samples only \({\lambda }_{4}\) showed a lwpp of more than 5 mT, while the linewidths of the other species matched the value of 2.5 mT quite well (Supporting Figure B10).

To best analyze the quality and robustness of both methods, the X-band cw-EPR spectrum and the pseudo-modulated and τ-averaged Q-band pulse EPR spectrum of sample Av1-Se2B-1 (Fig. 6) were analyzed by both methods, and the results were compared (Fig. 7). Because lwpp is a second independent parameter in the grid-of-errors approach, slightly better results are obtained in the simulation of X-band cw-EPR spectra than using the regularization (Fig. 7, upper panels). On the other hand, the lwpp parameter is slightly overestimated by the grid-of-errors method in Q-band (Fig. 7, lower panels), which lowers the quality of these results. Overall, spectral simulations obtained from either method are of excellent quality and show only minor deviations from the experimental data.

This detailed analysis hence demonstrates that regularization is a powerful and fast approach to simulate EPR spectra that are either dominated by only one statistically-distributed parameter (in this case, \(\lambda\)), or depend only on a second, non-dominant parameter (in this case, lwpp). To further improve the accuracy of the distribution P(\(\lambda\)), the samples could be measured in several frequency bands³⁷, and evaluated using a global regularization analogous to the analysis of DEER data sets³⁴. In summary, we believe that this method is more applicable to a high-spin EPR system than any simple strain models, since it provides faster, better, and model-free results for systems with many states and therefore many parameters.

--- Please insert Fig. 7 here ---

EPR analysis and assignment of Se-incorporated samples.

Pulsed- and cw-EPR experiments revealed that all species contain a total spin of 3/2, and all Se-species (\({\lambda }_{\text{1,3}-5}\)) relax faster than the FeMo resting state (\({\lambda }_{2}\)). 3P-ESEEM experiments of sample Av1-Se⁷⁷Se, which is labeled only at position 2B, confirms that Se is incorporated into the cofactor as its presence leads to additional hyperfine couplings. The same interpretation can be assumed for sample Av1-³³Se, although only spectroscopic, but no crystallographic confirmation is available for this sample²⁸. Spectral simulation revealed that the A_y value of the Se hyperfine coupling is 10.5 MHz, while the other two principal values A_x and A_z have to be treated with caution as their values only moderately influence the quality of the simulations. In addition to dead-time artifacts and cross suppression, the different hyperfine couplings of the individual spin species also impede unambiguous simulation results. The two S-labeled samples (Av1-S and Av1-³³S) were generated under turnover conditions in the presence of KSCN without N₂. Sample AvI-³³S exhibits additional hyperfine and quadrupole couplings with a strength of only a few MHz, which originate from one ³³S, and demonstrate that S exchange occurs even in the absence of N₂. We note that ENDOR experiments using uniformly ³³S-labeled nitrogenase have already been conducted. ³³S hyperfine couplings between − 10 MHz and − 16 MHz, including a quadrupole coupling of ~ 1 MHz, have been reported, but no specific S atom could be assigned¹⁹. In summary, additional hyperfine couplings in the 3P-ESEEEM spectra can be simulated by only one additional isotope (³³S or ⁷⁷Se).

All Se-incorporated samples contain four additional spin species (\({\lambda }_{\text{1,3}-5}\)), indicating that Se-exchange is possible under all experimental conditions studied (Table 1), most likely with yields above 90% at position 2B^28,30. Results from regularization (and from the grid-of-errors approach) show that regardless of the expected distribution of Se within the sulfur belt, the cw-EPR spectra always show similar rhombicity distributions and vary only in their probability intensities (Fig. 3D-H and Table 2). As crystallographic studies confirm different labeling pattern²⁸, it is possible that only the exchange at position 2B is detected spectroscopically and that additional Se exchange at positions 3A and 5A does not involve further changes in the electronic structure of the cluster. Note that Henthorn and colleagues carried out cw-EPR measurements with similarly prepared samples and detected no relevant changes in the EPR signals (Figure S2 in Ref. ³⁰). This does not contradict our results, as a closer look at their cw-EPR spectra reveals some additional low-intensity signals. Besides slightly different sample preparations, the reason could be the increased temperature of their measurements (10 K versus 5 K). Comparable cw-EPR measurements at 12 K support this interpretation: due to the short relaxation times of the FeMo cofactor, only a significantly broadened Se-pattern of low intensity can be detected (Supporting Figure B7).

To gain first insights into the nature of the four Se species, published EPR parameters of freeze-quenched reaction intermediates of the FeMo cofactor were extracted and compare with our values^14,22,23,49. Two identified \(S=3/2\) spin states ("1b” and "1c”)^22,49 have been previously assigned to hydride isomers of state E₂(2H)¹⁴. The effective \(g{\prime }\) factors of these species were extracted, and by using equations \(\lambda =\frac{2\left({\Delta }g\right)}{3{\left({\Delta }g\right)}^{2}-1}\), \({\Delta }g=\frac{{{g}_{\text{y}}^{{\prime }}}^{1/2}+{{g}_{\text{x}}^{{\prime }}}^{1/2}}{{{g}_{\text{y}}^{{\prime }}}^{1/2}-{{g}_{\text{x}}^{{\prime }}}^{1/2}}\) and assuming \({g}_{\text{x}}={g}_{\text{y}}\), the rhombicity values of these states may be calculated as \({\lambda }_{1\text{b}}\) ≈ 0.04 and \({\lambda }_{1\text{c}}\) ≈ 0.114, respectively. Other studies assigned a \(S=3/2\) spin state with a \(\lambda\) value of about 0.12 to the protonated resting state E₀(H⁺)^50,51, and a photoinduced state with a \(\lambda\) value of about 0.08 was very recently assigned to the (protonated) E₂ state⁵¹. Moreover, in freeze-quench experiments during turnover using an α-70Ile variant of Av1, a rhomboid signal (\(\lambda\) = 0.24) was assigned to state E₂⁵². This state can be excluded to be present in any of our samples as \(\lambda\) = 0.24 is well above the \(\lambda\) values observed in our spectra. The fact that a variant was used could explain the different \(\lambda\) values of the study and the one by Chica and coworkers⁵¹. An \(S\) = 1/2 signal with a \(g\)-factor of ≈ 2.00 that was assigned to state E₄^15,52 is again not observed in any of our Se-incorporated Av1 spectra.

These literature values and the λ values determined in this study are summarized in Table 3 and allow a first comparison: Species \({\lambda }_{1}\) has very similar values to the assigned state E₂(2H), species \({\lambda }_{3}\) to the assigned (hydrogenated) state E₂, and species \({\lambda }_{4}\) to the assigned states E₂(2H) or E₀(H⁺). Species \({\lambda }_{5}\) has never been observed in any other EPR experiment yet. Despite geometric distortions, \({\lambda }_{5}\) could represent another hydride isomer of E₂ or of any other higher state, as long as the total spin is S = 3/2.

The combination of the literature comparison, the analysis of the species distribution of sample Av1-S-remigration, and the results of the experiments with blue-light irradiation together support a more definite assignment of the different Se-species: Exchange of Se back to S, which was the rationale behind preparing sample Av1-S-remigration, does lead to a reduction of the states \({\lambda }_{1}\) and \({\lambda }_{4}\), but besides the ground state (\({\lambda }_{2}\)), state \({\lambda }_{3}\) persist even after prolonged reaction cycles. The light-irradiation experiments show very similar results: The probability of state \({\lambda }_{3}\) (and \({\lambda }_{2}\)) does not change in response to light. On the other hand, \({\lambda }_{1}\) and \({\lambda }_{4}\) respond reversibly to blue light, but whether a conversion of the hydrides upon light illumination really takes place, or only partial photolysis, needs to be clarified by further experiments.

Therefore, states \({\lambda }_{1}\) and \({\lambda }_{4}\), whose \(g\)-factors are very similar as those reported by¹⁴, and were referred to as states 1b and 1c, are most likely two different hydride isomers of state E₂. State \({\lambda }_{3}\), on the other hand, is irreversibly formed, representing a non-productive state that cannot be re-exchanged. It could either arise from a geometrical distortion due to the Se incorporation, or it could be a "stable" protonated E₀ state, which is irreversible due to the different pK_a value of the Se-FeMo cofactor (see below).

If the published assignments of the intermediate stages were not to be trusted, could in principle all Se species originate from geometric distortions? A number of findings speak against such an interpretation: First, it is highly unlikely that the \(g\)-values of unlabeled FeMo intermediates and of geometrically distorted Se-FeMo cofactors are very similar (Table 3). Second, if the individual Se-species would result from a simple geometric distortion of the FeMo cofactor by incorporation of the larger Se atom, either none or all of the Se-species would be re-exchanged by S in the Av1-S-remigration sample. Third, geometric distortions would lead to different \(\lambda\)-distributions of samples with Se-exchange at position 2B (samples Av1-Se2B-1 and Av1-Se2B-lowflux) compared to samples with an equal labeling of the sulfur belt (sample Av1-Se-C₂H₂), yet, our distributions do not show differences between the samples. In this context, it must be noted that the ground state values of the Se-FeMo and FeMo cofactors differ slightly (green dashed lines in Fig. 3), and hence potential differences in geometry may have a minor influence on the \(\lambda\) values.

A further important aspect in the discussion of the individual Se species is the decreased signal intensity of the Se-labeled samples compared to the S- or unlabeled samples: 40% of the FeMo cofactors are in an EPR silent state, which confirms that EPR inactive intermediate states of the FeMo cluster like E₁ or E₃ are also stabilized by the Se-incorporation method. As the S-labeled FeMo cofactors have the same cw-EPR intensity as the unlabeled cofactor, S-to-S exchange does not stabilize any intermediate states.

--- Please insert Table 3 here ---

Mechanistic insights. The question remains why intermediate states are stabilized by incorporation of Se. Basically, Se has a higher polarizability compared to S, and the Se-H group has a lower pK_a value compared to the S-H group⁵³, while serving as a structural surrogate for S in iron-sulfur clusters ⁵⁴. Moreover, calculations on Se (or S) metal model complexes discovered that the substitution of S with Se leads to a reduction of the ligand field strength and can additionally affect the energy of the electronic states ⁵⁵. These differences could lead to an equilibrium shift of the overall reaction upon Se substitution within the FeMo cofactor, and favor side reactions to early intermediate states (E₄, E₃, E₂, and E₁) accompanied by the release of H₂. No states higher than E₂ are observed, suggesting that the incorporation of Se into the FeMo cluster has to occur very early in the reaction scheme. The incorporation of Se into the 2B position of the FeMo cofactor could be accomplished via different reaction pathways^56,57. Our results support mechanisms that include protonation steps, as direct Se labeling would likely not result in as many different hydride isomers.

Earlier experiments with Se-modified samples showed that remigration of S results in a delayed enzyme activity²⁶, which lead to the assumption that the Se incorporated cofactor has a different activity and only regains its full enzymatic activity after S remigration. Our results demonstrate that Se incorporation leads to stabilization of different intermediate states containing different electronic structures. These differences could be due to changes in the effective oxidation states of the Fe atoms in the FeMo cofactor, whereby the total spin of S = 3/2 must be maintained. X-ray spectroscopy with a Se-labeled FeMo cofactor showed that position 2B and positions 3A/5A are electronically different³⁰. It was observed that the two iron atoms (Fe2/Fe6) that bind the Se at the 2B position show a "local oxidized character", whereas the iron nuclei which bind to the Se atoms at positions 3A/5A are rather reduced. It was also noted that both the incorporation of Se and hydrogen bonds affect the effective oxidation state and the electronic structure³⁰.

Can the additional protons of the E₂(2H) intermediate states be detected and characterized by EPR spectroscopy? Basically, additional protons show up in 3P-ESEEM spectra as additional signals around the proton Larmor frequency. Insets in Fig. 6B show these regions magnified for samples Av1-WT and Av1-Se2B-1. The proton hyperfine couplings of the Se-labeled sample (red) show a broadening compared to the unlabeled sample (black), and a weak splitting can be observed in the spectrum at the magnetic-field position 740 mT (see also Supporting Figure C4). Both of these indicate additional proton hyperfine couplings. ENDOR studies on the resting-state FeMo cofactor as well as on the CO-labeled cofactor have shown that the hyperfine couplings of the surrounding protons have only the strength of only a few MHz^19,58. It is therefore likely that any additional hyperfine couplings are only hardly visible in the 3P-ESEEM spectra due to fast relaxation times, low modulation depth, cross-suppression effects and are masked by the linewidth. It can still be concluded that the incorporation of Se leads to a broadened proton hyperfine signal pattern that most likely originate from additional protons attached to the Se-FeMo cofactor. Again, ENDOR spectroscopy at about 2 K could be helpful to further characterize these additional protons, in particular as the signals from species \({\lambda }_{1-5}\) are at least partially spectrally separated; a combination of blue-light illumination and orientation selection can further reduce the number of Se-species and enable unequivocal assignment.

In this study, various Se incorporation experiments into the catalytically active FeMo cofactor of a nitrogenase were investigated by EPR spectroscopy, as the property of such labels, e.g., their different reactivity, are far from being fully understood²⁸. Cw-EPR spectra of Se-incorporated samples showed complex signal patterns compared to unlabeled samples. Using Tikhonov regularization, applied for the first time in such a problem, it was possible to assign four different electronic states, each with a total spin of 3/2, differing only in their rhombicity. An independent second analysis using a grid-of-errors approach confirmed these results.

Using EPR parameters from already assigned intermediate states of the FeMo cofactor and irradiation experiments with blue light, one of the states could be assigned to the ground state (E₀) of the cofactor, and the other to (protonated) intermediate states (E₀(H⁺) and E₂(2H), see Table 3). Only one state \({\lambda }_{3}\), could potentially be stabilized by geometric distortions of the cofactor. By pulsed-EPR spectroscopy, small hyperfine couplings of ⁷⁷Se (and of ³³S) could be detected and spectrally simulated. These experiments confirmed the incorporation of at least one Se (or S) atom under turnover conditions.

The reason for the accumulation of the different reaction intermediates by Se incorporation presumably arises from the stabilization of these states due to the differences in polarizability and pK_a values between Se and S. As only "early" intermediates of the LT scheme were detected, the opening and incorporation of Se (and presumably also of other substrates) is very likely to proceed in the first steps of the reaction. It is also important to mention that 40% of the FeMo cofactors are in an EPR-silent state after Se-incorporation. Even if state E₁ is most probable of these states, higher odd states are also in principle possible.

The results presented here demonstrate that cw-EPR spectroscopy combined with Tikhonov regularization can analyze complex spectra with multiple species; this approach may be applied to other systems that contain paramagnetic transition metal centers. The result that under the selected experimental conditions a defined incorporation of Se or S takes place and reaction intermediates can be stabilized without effort offers great potential with respect to further investigations using different molecular spectroscopy methods.

continuous-wave (cw), Lowe and Thorneley (LT), MoFe-protein and Fe-protein from Azotobacter vinelandii (designated as Av1 and Av2), intrinsic peak-to-peak Lorentzian linewidth (lwpp), signal-to-noise ratio (S/N).

Author Contributions

L.H., K.P., T.S., O.E., S.W., D.R. and E.S. designed the research and conceived the experiments. K.P. and T.S. prepared all samples. L.H. conducted all EPR and ESEEM experiments. L.H. and E.S. analyzed and interpreted the spectroscopic data. L.H. wrote the simulation routines. The figures were generated by L.H. and E.S. The manuscript was written through contributions of all authors. All authors have reviewed and approved the manuscript.

Funding Sources

US National Institutes of Health grant GM045162 (D.C.R.)

Howard Hughes Medical Institute (D.C.R.)

DFG grant PP 1927 (project ID311061829), and RTG 1976 (project ID 235777276) (O.E.).

Acknowledgement

This work was supported by the Hans-Fischer-Gesellschaft (to E.S.). L.H., S.W. and E.S. thank the SIBW/DFG for financing EPR instrumentation that is operated within the MagRes Center of the University of Freiburg. We thank Jon Rittle for kindly preparing K⁷⁷SeCN, and Maximilian Mayländer and Johannes Ruhnke for help with the blue-light EPR experiments. Support from NIH Grant GM045162 and the Howard Hughes Medical Institute to D.C.R. is gratefully acknowledged. O.E. acknowledges support from the German Research Foundation (PP 1927, project ID 311061829, and RTG 1976, project ID 235777276).

Boyd, E. S. & Peters, J. W. New insights into the evolutionary history of biological nitrogen fixation. Frontiers in Microbiology 4, 1-12 (2013).
Hoffman, B. M., Lukoyanov, D., Yang, Z.-Y., Dean, D. R. & Seefeldt, L. C. Mechanism of nitrogen fixation by nitrogenase: The next stage. Chemical Reviews 114, 4041-4062 (2014).
Einsle, O. & Rees, D. C. Structural enzymology of nitrogenase enzymes. Chemical Reviews 120, 4969-5004 (2020).
Kim, J. S. & Rees, D. C. Crystallographic Structure and Functional Implications of the Nitrogenase Molybdenum Iron Protein from Azotobacter-Vinelandii. Nature 360, 553-560 (1992).
Einsle, O. et al. Nitrogenase MoFe-protein at 1.16 Å resolution: A central ligand in the FeMo-cofactor. Science 297, 1696-1700 (2002).
Spatzal, T. et al. Evidence for interstitial carbon in nitrogenase FeMo cofactor. Science 334, 940 (2011).
Seefeldt, L. C., Hoffman, B. M. & Dean, D. R. Electron transfer in nitrogenase catalysis. Current Opinion in Chemical Biology 16, 19-25 (2012).
Rutledge, H. L. & Tezcan, F. A. Electron transfer in nitrogenase. Chemical Reviews 120, 5158-5193 (2020).
Hageman, R. V. & Burris, R. H. Kinetic studies on electron transfer and interaction between nitrogenase components from Azotobacter vinelandii. Biochemistry 17, 4117-4124 (1978).
Milton, R. D. & Minteer, S. D. Nitrogenase bioelectrochemistry for synthesis applications. Accounts of Chemical Research 52, 3351-3360 (2019).
Lowe, D. J. & Thorneley, R. N. F. The mechanism of Klebsiella pneumoniae nitrogenase action. The determination of rate constants required for the simulation of the kinetics of N₂ reduction and H₂evolution. Biochemical Journal 224, 895-901 (1984).
Hoffman, B. M., Lukoyanov, D., Dean, D. R. & Seefeldt, L. C. Nitrogenase: A draft mechanism. Accounts of Chemical Research 46, 587-595 (2013).
Hoeke, V. et al. High-resolution ENDOR spectroscopy combined with quantum chemical calculations reveals the structure of nitrogenase Janus intermediate E₄(4H). Journal of the American Chemical Society 141, 11984-11996 (2019).
Lukoyanov, D. A. et al. Hydride conformers of the nitrogenase FeMo-cofactor two-electron reduced state E2(2H), assigned using cryogenic intra electron paramagnetic resonance cavity photolysis. Inorganic Chemistry 57, 6847-6852 (2018).
Lukoyanov, D. et al. Reductive elimination of H₂ activates nitrogenase to reduce the N≡N triple bond: Characterization of the E₄(4H) Janus intermediate in wild-type enzyme. Journal of the American Chemical Society 138, 10674-10683 (2016).
Lukoyanov, D. et al. ENDOR/HYSCORE studies of the common intermediate trapped during nitrogenase reduction of N₂H₂, CH₃N₂H, and N₂H₄ support and alternating reaction pathway for N₂ reduction. Journal of the American Chemical Society 133, 11655-11664 (2011).
Cao, L. & Ryde, U. Putative reaction mechanism of nitrogenase after dissociation of a sulfide ligand. Journal of Catalysis 391, 247-259 (2020).
Rawlings, J. et al. Novel metal cluster in the iron-molybdenum cofactor of nitrogenase. Spectroscopic evidence. Journal of Biological Chemistry 253, 1001-1004 (1978).
Venters, R. A. et al. ENDOR of the resting state of nitrogenase molybdenum–iron proteins of Azotobacter vinelandii, Klebsiella pneumoniae, and Clostridium pasteurianum: ¹H, ⁵⁷Fe, ⁹⁵Mo, and ³³S studies. Journal of the American Chemical Society 108, 3487-3498 (1986).
Goldberg, D. P. et al. EPR spectra from “EPR-silent” species: High-field EPR spectroscopy of manganese(III) porphyrins. Journal of the American Chemical Society 119, 8722-8723 (1997).
Telser, J. EPR interactions – Zero-field splittings. eMagRes 6, 207-233 (2017).
Fisher, K. et al. Conformations generated during turnover of the Azotobacter vinelandii nitrogenase MoFe protein and their relationship to physiological function. Journal of Inorganic Biochemistry 101, 1649-1656 (2007).
Fisher, K., Newton, W. E. & Lowe, D. J. Electron paramagnetic resonance analysis of different Azotobacter vinelandii nitrogenase MoFe-protein conformations generated during enzyme turnover: Evidence for S=3/2 spin states from reduced MoFe-protein intermediates. Biochemistry 40, 3333-3339 (2001).
Shaw, S. et al. Nitrite and hydroxylamine as nitrogenase substrates: Mechanistic implications for the pathway of N₂ reduction. Journal of the American Chemical Society 136, 12776-12783 (2014).
Benton, P. M. C. et al. Localization of a substrate binding site on the FeMo-cofactor in nitrogenase: Trapping propargyl alcohol with an α-70-substituted MoFe protein. Biochemistry 42, 9102-9109 (2003).
Spatzal, T., Perez, K. A., Einsle, O., Howard, J. B. & Rees, D. C. Ligand binding to the FeMo-cofactor: structures of CO-bound and reactivated nitrogenase. Science 345, 1620-1623 (2014).
Rohde, M., Grunau, K. & Einsle, O. CO binding to the FeV cofactor of CO-reducing vanadium nitrogenase at atomic resolution. Angewandte Chemie International Edition 59, 23626-23630 (2020).
Spatzal, T., Perez, K. A., Howard, J. B. & Rees, D. C. Catalysis-dependent selenium incorporation and migration in the nitrogenase active site iron-molybdenum cofactor. eLife 4, e11620 (2015).
Sippel, D. et al. A bound reaction intermediate sheds light on the mechanism of nitrogenase. Science 359, 1484-1489 (2018).
Henthorn, J. T. et al. Localized electronic structure of nitrogenase FeMoco revealed by selenium K-edge high resolution X-ray absorption spectroscopy. Journal of the American Chemical Society 141, 13676-13688 (2019).
Spatzal, T. et al. Nitrogenase FeMoco investigated by spatially resolved anomalous dispersion refinement. Nature Communications 7, 10902 (2016).
Chiang, Y.-W., Borbat, P. P. & Freed, J. H. The determination of pair distance distributions by pulsed ESR using Tikhonov regularization. Journal of Magnetic Resonance 172, 279-295 (2005).
Jeschke, G. et al. DeerAnalysis2006 – a comprehensive software package for analyzing pulsed ELDOR data. Applied Magnetic Resonance 30, 473-498 (2006).
Rein, S., Lewe, P., Andrade, S. L., Kacprzak, S. & Weber, S. Global analysis of complex PELDOR time traces. Journal of Magnetic Resonance 295, 17-26 (2018).
Edwards, T. H. & Stoll, S. Optimal Tikhonov regularization for DEER spectroscopy. Journal of Magnetic Resonance 288, 58-68 (2018).
Ibáñez, L. F. & Jeschke, G. General regularization framework for DEER spectroscopy. Journal of Magnetic Resonance 300, 28-40 (2019).
Azarkh, M. & Groenen, E. J. J. Simulation of multi-frequency EPR spectra for a distribution of the zero-field splitting. Journal of Magnetic Resonance 255, 106-113 (2015).
Stoll, S. & Schweiger, A. EasySpin, a comprehensive software package for spectral simulation and analysis in EPR. Journal of Magnetic Resonance 178, 42-55 (2006).
Stoll, S. & Britt, R. D. General and efficient simulation of pulse EPR spectra. Physical Chemistry Chemical Physics 11, 6614-6625 (2009).
Lee, H.-I., Doan, P. E. & Hoffman, B. M. General analysis of ¹⁴N (I = 1) electron spin echo envelope modulation. Journal of Magnetic Resonance 140, 91-107 (1999).
George, G. N., Prince, R. C. & Bare, R. E. Electron paramagnetic resonance spectroscopy of the iron−molybdenum cofactor of Clostridium pasteurianum nitrogenase. Inorganic Chemistry 35, 434-438 (1996).
Hagen, W. R. Wide zero field interaction distributions in the high-spin EPR of metalloproteins. Molecular Physics 105, 2031-2039 (2007).
Hagen, W. R., Hearshen, D. O., Harding, L. J. & Dunham, W. R. Quantitative numerical analysis of g strain in the EPR of distributed systems and its importance for multicenter metalloproteins. Journal of Magnetic Resonance 61, 233-244 (1985).
Hearshen, D. O. et al. An analysis of g strain in the EPR of two [2Fe-2S] ferredoxins. Evidence for a protein rigidity model. Journal of Magnetic Resonance 69, 440-459 (1986).
Froncisz, W. & Hyde, J. S. Broadening by strains of lines in the g‐parallel region of Cu²⁺ EPR spectra. The Journal of Chemical Physics 73, 3123-3131 (1980).
Lee, H.-I. et al. The interstitial atom of the nitrogenase FeMo-cofactor: ENDOR and ESEEM show it is not an exchangeable nitrogen. Journal of the American Chemical Society 125, 5604-5605 (2003).
Azarkh, M., Gast, P., Mason, A. B., Groenen, E. J. J. & Mathies, G. Analysis of the EPR spectra of transferrin: The importance of a zero-field-splitting distribution and 4th-order terms. Physical Chemistry Chemical Physics 21, 16937-16948 (2019).
Nehrkorn, J., Telser, J., Holldack, K., Stoll, S. & Schnegg, A. Simulating frequency-domain electron paramagnetic resonance: bridging the gap between experiment and magnetic parameters for high-spin transition-metal ion complexes. Journal of Physical Chemistry B 119, 13816-13824 (2015).
Lukoyanov, D. et al. A confirmation of the quench-cryoannealing relaxation protocol for identifying reduction states of freeze-trapped nitrogenase intermediates. Inorganic Chemistry 53, 3688-3693 (2014).
Morrison, C. N., Spatzal, T. & Rees, D. C. Reversible protonated resting state of the nitrogenase active site. Journal of the American Chemical Society 139, 10856-10862 (2017).
Chica, B. et al. Defining intermediates of nitrogenase MoFe protein during N₂ reduction under photochemical electron delivery from CdS quantum dots. Journal of the American Chemical Society 142, 14324-14330 (2020).
Lukoyanov, D., Barney, B. M., Dean, D. R., Seefeldt, L. C. & Hoffman, B. M. Connecting nitrogenase intermediates with the kinetic scheme for N₂ reduction by a relaxation protocol and identification of the N₂ binding state. Proceedings of the National Academy of Sciences of the United States of America 104, 1451-1455 (2007).
Wessjohann, L. A., Schneider, A., Abbas, M. & Brandt, W. Selenium in chemistry and biochemistry in comparison to sulfur. Biological Chemistry 388, 997-1006 (2007).
Zheng, B., Chen, X.-D., Zheng, S.-L. & Holm, R. H. Selenium as a structural surrogate of Sulfur: Template-assisted assembly of five types of Tungsten–Iron–Sulfur/Selenium clusters and the structural fate of chalcogenide reactants. Journal of the American Chemical Society 134, 6479-6490 (2012).
Spiller, N., Chilkuri, V. G., DeBeer, S. & Neese, F. Sulfur vs. selenium as bridging ligand in di-iron complexes: A theoretical analysis. European Journal of Inorganic Chemistry 2020, 1525-1538 (2020).
Arias, R. J. Examination of selenium incorporation and product formation in the nitrogenase FeMo-cofactor. Ph.D thesis, Califorina Institute of Technology Pasadena California (USA) (2018).
Dance, I. Mechanisms of the S/CO/Se interchange reactions at FeMo-co, the active site cluster of nitrogenase. Dalton Transactions 45, 14285-14300 (2016).
Lee, H.-I., Cameron, L. M., Hales, B. J. & Hoffman, B. M. CO binding to the FeMo Cofactor of CO-inhibited nitrogenase: ¹³CO and ¹H Q-band ENDOR investigation. Journal of the American Chemical Society 119, 10121-10126 (1997).

Table 1

MoFe-protein samples and modifications used in this study. All turnover assays are performed without nitrogen.
Sample	Abbreviation	Sample condition	Labeling (based on ^28,30)	Concentration
Azotobacter vinelandii MoFe protein (Av1) wild type – resting state	Av1-WT	100% Av1	--	≈ 48 mg/mL
Av1 wild type (turnover)	Av1-Se2B-1	10 mM KSeCN was added to the turnover assay (Av2/Av1 ratio = 2).	Se (position 2B)	≈ 73 mg/mL
Av1 wild type (turnover)	Av1-⁷⁷Se2B	10 mM K⁷⁷SeCN was added to the turnover assay (Av2/Av1 ratio = 2).	⁷⁷Se (position 2B)	≈ 74 mg/mL
Av1 wild type (turnover)	Av1-Se-low	0.25 mM KSeCN was added to the turnover assay (Av2/Av1 ratio = 2).	Se (predominantly at position 2B)	≈ 70 mg/mL
Av1 wild type (turnover)	Av1-³³S	10 mM K³³SCN was added to the turnover assay (Av2/Av1 ratio = 2).	³³S (position 2B)	≈ 40 mg/mL
Av1 wild type (turnover)	Av1-Se-C₂H₂	sample Av1-Se2B-1 was used for a second turnover assay, which was quenched with 10 mM KSeCN after a reaction time of 5 min.	Se (positions 2B, 3A and 5A)	≈ 70 mg/mL
Av1 wild type (turnover)	Av1-Se2B-lowflux	15 mM KSeCN was added to the turnover assay (Av2/Av1 ratio = 1.5).	Se (position 2B)	≈ 57 mg/mL
Av1 wild type (turnover)	Av1-S	22.5 mM KSCN was added to the turnover assay (Av2/Av1 ratio = 1.5).	S (position 2B)	≈ 40 mg/mL
Av1 wild type (prolonged turnover)	Av1-S-remigration	sample Av1-Se2B-lowflux was used. A second turnover assay was started with an Av2/Av1 ratio = 4. The assay proceeded for ≈ 1 h.	Se is expected to be replaced by S again.	≈ 46 mg/mL

Table 2

Summary of rhombicity parameters extracted from Fig. 3. \({\lambda }\) values were extracted manually by determining the local maxima of the distribution, values with a probability height of less than 10% were ignored from further discussion. The respective populations of all Se-incorporated samples were fitted using a multi-Gaussian function (Supporting Table B11), and all fractions were determined. The respective largest fraction is depicted in bold. Effective \({g}\)-factors were calculated by using Eq. 1-SI.
		\({{\lambda }}_{1}\) (fraction %)	\({{\lambda }}_{2}\) (fraction %)	\({{\lambda }}_{3}\) (fraction %)	\({{\lambda }}_{4}\) (fraction %)	\({{\lambda }}_{5}\) (fraction %)
Avl-WT			0.055
Av1-S			0.053
Av1-³³S			0.053
Av1-Se2B-1		0.033 (11.5%)	0.056 (17.2%)	0.080 (23.3%)	0.113 (35.9%)	0.166–0.190 (12.1%)
Av1-Se2B-lowflux		0.033 (14.4%)	0.058 (34.3%)	0.086 (18.4%)	0.118 (30%)	≈ 0.189 (2.9%)
Av1-⁷⁷Se2B		0.032 (14.1%)	0.056 (20.2%)	0.080 (10.6%)	0.115 (37.6%)	0.163–0.220 (4.3%)
Av1-Se-low		0.033 (24.2%)	0.058 (8.4%)	0.078 (22%)	0.120 (38.9%)	0.171–0.230 (4.3%)
Av1-Se-C₂H₂		0.034 (19.2%)	0.058 (20.8%)	0.084 (10.6%)	0.120 (33.4%)	0.175–0.217 (6.5%)
Av1-S-remigration		0.033 (4.1%)	0.057 (30.5%)	0.086 (39.4%)	0.112 (16.1%)	0.180 (9.9%)
Average value		0.033	0.056	0.082	0.116	≈ 0.19
Effective g-values*	\({{g}^{{\prime }}}_{\text{x}}^{1/2}\) \({g{\prime }}_{\text{y}}^{1/2}\) \({g{\prime }}_{\text{z}}^{1/2}\)	3.80 4.20 2.03	3.66 4.33 2.02	3.50 4.48 2.02	3.30 4.68 2.00	≈ 2.92 ≈ 4.98 ≈ 1.82
* Calculated from Supporting Information Eq. 1 and \({g}_{x}={g}_{y}=2.00\)and\({g}_{z}=2.03\)

There is NO Competing Interest.

NatCommNitroSeSIfinalV2.docx
TOC.png
TOC figure

Download PDF

Journal Publication

published 13 May, 2024

Read the published version in Nature Communications →

Version 1

posted

You are reading this latest preprint version

Analysis of early intermediate states of the nitrogenase reaction by Se incorporation and regularization of EPR spectra

Status:

Journal Publication

Version 1

Abstract

Figures

Introduction

Experimental Section

Results

Discussion

Summary and Outlook

Abbreviations

Declarations

References

Tables

Additional Declarations

Supplementary Files

Status:

Journal Publication

Version 1