Modeling

Overview

Mathematical modeling utilizes concrete numbers, estimations, simulations, assumptions, descriptive equations and statistics to determine the outcome or predictable results of a particular system of operation. Mathematical modeling takes real life problems and examples and converts them into a series of descriptive equations to solve for. Through the process of modeling the equations and defining various variables that we call parameters, we are able to conduct analysis regarding the relationship—provided one exsists—between several variables and their individual and independent effects on the system. We will be able to see the results of the presence or the absence of the particular variable to deem its importance to the model that we’ve constructed to simulate a real life situation. We can use it to validate our hypotheses or the success of our lab experiments. Modeling an experiment before conducting it can help researchers optimize the concentration of reactants they may need to achieve a result they are looking for.



For our team, we used MATLAB to model differential equations including DNA transcription, translation, protein folding, dimerization, tetramerization, and the rate of protein association/dissociation to the promoter as well as protein degradation.







Genetic Circuit

Oscillatory expression of a target protein can be achieved by integrating its correlated gene into a regulatory network involving positive and negative feedback. A hybrid promoter that allows for simultaneous transcription of activator, araC, and repressor, lacI, components creates a dual feedback circuit. Production of the AraC protein results in increased promoter activity causing increased expression of LacI. The circuit then undergoes negative feedback, decreasing promoter activity. We created a computational model including the rates of transcription, translation, multimerization, protein folding, and DNA looping. This provides predictions of the oscillatory expression of a target gene. This system can be altered to regulate the expression of varied target genes, making it a useful tool in synthetic biology.


We varied IPTG concentration and arabinose percentage to achieve the robust oscillations of our protein of interest, which we chose AraC dimer molecules for. In addition, we wrote differential equations to demonstrate the impact of incorporating our nitrogenase cluster into the genetic circuit, and to verify the robustness of the oscillations.

We also modeled and compared two versions of the code. The first contained three hyrbid promoters, which control the transcription of araC, lacI and gfp. The second contained 4 hyrbid promoters, controlling the transcription of araC, lacI, gfp and the nifH gene. We found that proper nifH transcription was critical to the assembly of a fuctional nitrogenase enzyme [1]. We adjusted IPTG and Arabinose concentrations to determine the initial concentration we should use to induce our cells to generate robust oscillations with varying oscillation periods

We worked on the MATLAB code alongside our wet lab experiments by:

  • Using our model to determine the IPTG concentration and arabinose percentage that needs to be included in the system to ensure robust oscillation

  • Determining whether a three or four promoter system supports robust oscillations and how the increase in promoters affects the genetic circuit

  • Predicting levels of protein expression at different moments of time

  • Utilizing wet lab experiments for verification of oscillatory behavior

  • Confirming period times match the ones our model predicts for different initial concentrations of IPTG and arabinose





Details of Our Model

Our model improves on attempts of modeling Hasty’s 2008 Oscillator. Our model is the only iGEM team model that utilizes mass ODEs to predict oscillatory behavior of this system. We chose to use mass ODEs, where we define each step in the process, because this mechanistic approach is simple to understand and replicable for other teams. Past teams, like the 2015 HZAU China team, utilized delayed differential equations to model the genetic oscillator to correct an incomplete model. These flaws likely originated from their incorrectly defined x term, which represents the number of ssrA tags in the system. They defined x as a constant throughout the entire code, which affects the ability of the genetic circuit to oscillate. In contrast, we defined x as a variable that changed accordingly to the concentration of each protein with an ssrA tag.

Using our model, we were able to understand how the oscillations were impacted by the concentrations of IPTG and arabinose that we used to induce our system. Based on our model, we identified the combinations of IPTG and arabinose concentration that produced robust oscillations. After narrowing down the initial conditions that gave robust oscillations, we created a table of the average periods of each oscillation, and the average amplitude of the corresponding oscillations. Using this information, we tuned our genetic oscillator in experiments to produce short oscillations, long oscillations, low amplitudes, and high amplitudes.

To verify the accuracy of our model, we performed confocal microscopy to visualize GFP oscillations. We induced our cells at 0.7% arabinose and 2 mM IPTG and compared the experimental results to our oscillations in the model. The cells were nearly instantly visible with green fluorescence after induction, as predicted by our model. We were able to obtain cells that had a peak oscillation after around an hour, which matches what is predicted by our model after the initial oscillation.

Our modeling approach provides a good example for others because it utilizes simple nomenclature and a logical flow to predict protein expression. By modeling each part of the translation process using ODEs, we separated our code into transcription, translation, protein folding, dimerization/tetramerization, and protein association sections. This is an intuitive format that is easy to follow along with our flowchart of the mechanism of the genetic oscillator. Additionally, our code is thoroughly annotated to avoid any confusion in nomenclature or any other parts of the code.






Comparison

Arabinose (%) and IPTG (mM) were found to be important in controlling period, amplitude, and robustness of oscillations [2]. Using the equations below, period and amplitudes for stable oscillations were determined for different arabinose (%) and IPTG (mM) levels over an interval of 600 minutes.

$$findpeaks(specs(:,9),MinPeakHeight=100);$$ $$meanperiod = mean(diff(t1(locs)))$$


In the following tables, red indicates dying oscillations (not robust), yellow indicates lowest value among robust oscillations, and green indicates highest value among robust oscillations.


Period (min) among different values of arabinose (%) and IPTG (mM) in 3 promoter model.


Controlling for IPTG concentration, arabinose percentage of 0.8 produced the longest period on average among robust oscillations and arabinose percentage of 0.1 produced the shortest period on average among robust oscillations.


Controlling for arabinose, IPTG concentration of 1.2 mM produced the longest period on average among robust oscillations and IPTG concentration of 0 mM produced the shortest period on average among robust oscillations.


Arabinose percentage of 0.1 and IPTG concentration of 0 mM demonstrated the shortest period among robust oscillations.


Arabinose percentage of 0.9 and IPTG concentration of 1.8 demonstrated the longest period among robust oscillations.



Period (min) among different values of arabinose (%) and IPTG (mM) in 4 promoter model.


As arabinose increased, the number of possible IPTG concentrations that would result in a robust oscillation decreased.


Controlling for IPTG concentration, arabinose percentage of 0.4 produced the longest period on average among robust oscillations and arabinose percentage of 0.1 produced the shortest period on average among robust oscillations.


Controlling for arabinose, IPTG concentration of 1.2 mM produced the longest period on average among robust oscillations and IPTG concentration of 0 mM produced the shortest period on average among robust oscillations.


Arabinose percentage of 0.1 and IPTG concentration of 0 mM demonstrated the shortest period among robust oscillations. Arabinose percentage of 0.4 and IPTG concentration of 1.6 mM demonstrated the longest period among robust oscillations.



These tables were used to generate graphs that better can demonstrate the differences between the 3 promoter model and 4 promoter model. Only values that produced robust oscillations were included in the graphs.



3 Promoter Model Periods (Oscillator without Load). Period (min) vs. IPTG (mM), Arabinose (%).



4 Promoter Model Periods (Oscillator with Load). Period (min) vs. IPTG (mM), Arabinose (%).


The 3 promoter model demonstrates more tunability such that it shows more possibilities of robust oscillations under the combinations of arabinose and IPTG tested.


3 Promoter Model Output (Oscillator without Load). Arabinose 0.4%. IPTG 1.8 mM.


4 Promoter Model Output (Oscillator with Load). Arabinose 0.4%. IPTG 1.8 mM.


3-Promoter and 4-Promoter Model— Nomenclature

mA araC mRNA molecules
mR lacI mRNA molecules
mG GFP mRNA molecules
mNifH nifH mRNA molecules
AraCuf unfolded AraC monomers
LacIuf unfolded LacI monomers
GFPuf unfolded GFP monomers
AraC1 AraC monomers (folded)
AraC2 AraC dimers
LacI1 LacI monomers (folded)
LacI2 LacI dimers
LacI4 LacI tetramers
GFP GFP monomers (folded)
NifH1 NifH monomers
NifH2 NifH dimers
p00A araC promoters with no proteins bound
p01A araC promoters with 1 LacI tetramer bound
p02A araC promoters with 2 LacI tetramers bound
p10A araC promoters with 1 AraC dimer bound
p11A araC promoters with 1 AraC dimer bound and 1 LacI tetramer bound
p12A araC promoters with 1 AraC dimer bound and 2 LacI tetramers bound
pL0A araC promoters that are looped with no LacI tetramers bound
pL1A araC promoters that are looped with 1 LacI tetramer bound
pL2A araC promoters that are looped with 2 LacI tetramers bound
p00L lacI promoters with no proteins bound
p01L lacI promoters with 1 LacI tetramer bound
p02L lacI promoters with 2 LacI tetramers bound
p10L lacI promoters with 1 AraC dimer bound
p11L lacI promoters with 1 AraC dimer bound and 1 LacI tetramer bound
p12L lacI promoters with 1 AraC dimer bound and 2 LacI tetramers bound
pL0L lacI promoters that are looped with no proteins bound
pL1L lacI promoters that are looped with 1 LacI tetramer bound
pL2L lacI promoters that are looped with 2 LacI tetramers bound
p00G gfp promoters with no proteins bound
p01G gfp promoters with 1 LacI tetramer bound
p02G gfp promoters with 2 LacI tetramers bound
p10G gfp promoters with 1 AraC dimer bound
p11G gfp promoters with 1 AraC dimer bound and 1 LacI tetramer bound
p12G gfp promoters with 1 AraC dimer bound and 2 LacI tetramers bound
pL0G gfp promoters that are looped with no proteins bound
pL1G gfp promoters that are looped with 1 LacI tetramer bound
pL2G gfp promoters that are looped with 2 LacI tetramers bound
p00N gfp promoters with no proteins bound
p01N gfp promoters with 1 LacI tetramer bound
p02N gfp promoters with 2 LacI tetramers bound
p10N gfp promoters with 1 AraC dimer bound
p11N gfp promoters with 1 AraC dimer bound and 1 LacI tetramer bound
p12N gfp promoters with 1 AraC dimer bound and 2 LacI tetramers bound
pL0N gfp promoters that are looped with no proteins bound
pL1N gfp promoters that are looped with 1 LacI tetramer bound
pL2N gfp promoters that are looped with 2 LacI tetramers bound

3-Promoter and 4-Promoter Model— Parameters

Arabinose 0.4 Arabinose concentration in %
IPTG 1.6 IPTG concentration in mM
bA 0.36 araC transcription rate in min-1
bR 0.36 lacI transcription rate in min-1
bG 0.36 GFP transcription rate in min-1
bN 0.36 nif transcription rate in min-1
tA 90 AraC translation rate in min-1
tR 90 LacI translation rate in min-1
tG 90 GFP translation rate in min-1
tNifH 72.1176 NifH translation rate in min-1. Based on proportionality of GFP and NifH base pairs
kA 1.8 * (((arabinose^2) /((arabinose^2) + (2.5^2))) * (1/ (1 + (IPTG/1.8)^2))) Binding rate of AraC to the promoter
k-A 1.8 Dissociation rate of AraC from the promoter
kR 1.8 * ((0.2-0.01) * (1/ (1 + (IPTG/0.035)^2)) + 0.01) Binding rate of LacI to the promoter
k-R 1.8 Dissociation rate of LacI from the promoter
kFA 0.9 Folding rate of AraC monomer
kFR 0.9 Folding rate of LacI monomer
kFG 0.9 Folding rate of GFP monomer
kDA 0.018 Dimerization rate of AraC
k-DA 0.00018 Undimerization rate of AraC
kDR 0.018 Dimerization rate of LacI
k-DR 0.00018 Undimerization rate of LacI
kDNifH 0.018 Dimerization rate of NifH in min-1
k-DNifH 0.00018 Undimerization rate of NifH in min-1
kT 0.018 Tetramerization rate of LacI
k-T 0.00018 Untetramerization rate of LacI
kL 0.36 Rate of promoter DNA looping
kUL 0.18 Rate of promoter DNA unlooping
dA 0.54 Degradation rate of araC mRNA in min-1
dR 0.54 Degradation rate of lacI mRNA in min-1
dG 0.54 Degradation rate of GFP mRNA in min-1
dN 0.54 Degradation rate of Nif mRNA in min-1
x AraCuf + LacIuf + GFPuf + AraC1 + LacI1 + GFP + NifH1 + 2*AraC2 + 2*LacI2 + 4*LacI4 + 2*NifH2 + 4*(p01A + p01L + p01G + p01N) + 8*(p02A+p02L+p02G+p02N) + 2*(p10A + p10L + p10G + p10N) + 6*(p11A + p11L + p11G + p11N) + 10*(p12A + p12L + p12G + p12N) + 8*(pL2A + pL2L + pL2G + pL2N) + 4*(pL1A + pL1L + pL1G + pL1N) Number of ssRA tags
f_x 1080/(0.1 + x) Degradation rate for all the proteins

The equations related to nitrogenase are only present in the four promoter model. MATLAB code for the system is embedded below:

Change in mRNA over time:
$$ \frac{dMA}{dt} = b_A*p00A + 20*b_A*p10A - d_A*mA $$ $$ \frac{dMR}{dt} = b_R*p00L + 20*b_R*p10L - d_R*mR $$ $$ \frac{dMG}{dt} = b_G*p00G + 20*b_G*p10G - d_G*mG $$ $$ \frac{dMNifH}{dt} = b_N*p00N + 20*b_N*p10N - d_N*mNifH $$

Transcription is 20 times faster when there is an AraC dimer bound to the promoter. In the equation pxyZ, x represents the number of AraC dimers bound to the promoter, y represents the number of LacI tetramers bound to the promoter, and Z represents which promoter. x can be 0 or 1, y can be 0, 1, or 2, and Z can be A, R, G or N. LacI prevents the p00A and p10A states from being present in the system.


Change in unfolded protein over time:
$$ \frac{dAraC_{uf}}{dt} = t_A*mA - k_{FA}*AraC_{uf} - 2.5*f_x*AraC_{uf} $$ $$ \frac{dLacI_{uf}}{dt} = t_R*mR - k_{FR}*LacI_{uf} - f_x*LacI_{uf} $$ $$ \frac{dGFP_{uf}}{dt} = t_G*mG - k_{FG}*GFP_{uf} - 2.5*f_x*GFP_{uf} $$


Change in AraC proteins over time:
$$ \frac{dAraC_1}{dt} = k_{FA}*AraC_{uf} - 2*k_{d_A}*(AraC_1^2) + 2*k_{-d_A} $$ $$ \ *AraC_2 - 2.5*f_x*AraC_1 $$ $$ \frac{dAraC_2}{dt} = k_{d_A}*(AraC_1^2) - k_{-d_A}*AraC_2 - k_A*p00A*AraC_2 $$ $$ \ + k_{-A}*p10A - k_A*p01A*AraC_2 + k_{-A}*p11A - k_A*p02A*AraC_2 $$ $$ \ + k_{-A}*p12A + k_L*p12A - k_A*p00L*AraC_2 + k_{-A}*p10L $$ $$ \ - k_A*p01L*AraC_2 + k_{-A}*p11L - k_A*p02L*AraC_2 + k_{-A}*p12L $$ $$ \ + k_L*p12L - 2.5*f_x*AraC_2 - k_A*AraC_2*p00G + k_{-A}*p10G $$ $$ \ - k_A*AraC_2*p01G + k_{-A}*p11G - k_A*AraC_2*p02G + k_{-A}*p12G $$ $$ \ + k_L*p12G $$


Change in LacI proteins over time:
$$ \frac{dLacI_1}{dt} = k_{FR}*LacI_{uf} - 2*k_{d_R}*(LacI_1^2) $$ $$ \ + 2*k_{-d_R}*LacI_2 - f_x*LacI_1 $$ $$ \frac{dLacI_2}{dt} = k_{d_R}*(LacI_1^2) - k_{-d_R}*LacI_2 $$ $$ \ - 2*k_T*(LacI_2^2) + 2*k_{-T}*LacI_4 - f_x*LacI_2 $$ $$ \frac{dLacI_4}{dt} = k_T*(LacI_2^2) - k_{-T}*LacI_4 $$ $$ \ - 2*k_R*p00A*LacI_4 + k_{-R}*p01A - 2*k_R*p10A*LacI_4 $$ $$ \ + k_{-R}*p11A - k_R*p01A*LacI_4 + 2*k_{-R}*p02A $$ $$ \ - k_R*p11A*LacI_4 + 2*k_{-R}*p12A - 2*k_R*p00L*LacI_4 $$ $$ \ + k_{-R}*p01L - 2*k_R*p10L*LacI_4 + k_{-R}*p11L $$ $$ \ - k_R*p01L*LacI_4 + 2*k_{-R}*p02L - k_R*p11L*LacI_4 $$ $$ \ + 2*k_{-R}*p12L - f_x*LacI_4 - k_R*LacI_4*p00G + k_{-R}*p01G $$ $$ \ - k_R*LacI_4*p01G + k_{-R}*p02G - k_R*LacI_4*p10G $$ $$ \ + k_{-R}*p11G - k_R*LacI_4*p11G + k_{-R}*p12G $$


Change in GFP protein over time:
$$ \frac{dGFP}{dt} = k_{FG}*GFP_{uf} - 2.5*f_x*GFP $$


Change in NifH proteins over time:
$$ \frac{dNifH_1}{dt} = t_{NifH}*mNifH - 2*k_{d_NifH}*(NifH_1^2) $$ $$ \ + k_{-d_NifH}*NifH_2 - 2.5*f_x*NifH_1 $$ $$ \frac{dNifH_2}{dt} = 2*k_{d_NifH}*(NifH_1^2) $$ $$ \ - k_{-d_NifH}*NifH_2 - 2.5*f_x*NifH_2 $$


Change in araC promoter states over time:
$$ \frac{dP00A}{dt} = -k_A*p00A*AraC_2 + k_{-A}*p10A - 2*k_R*p00A*LacI_4 $$ $$ \ + k_{-R}*p01A + k_UL*pL0A + f_x*p10A + f_x*p01A $$ $$ \frac{dP01A}{dt} = -k_A*p01A*AraC_2 + k_{-A}*p11A + 2*k_R*p00A*LacI_4 $$ $$ \ - k_{-R}*p01A - k_R*p01A*LacI_4 + 2*k_{-R}*p02A + f_x*p11A $$ $$ \ - f_x*p01A + 2*f_x*p02A $$ $$ \frac{dP02A}{dt} = -k_A*p02A*AraC_2 + k_{-A}*p12A + k_R*p01A*LacI_4 $$ $$ \ - 2*k_{-R}*p02A - k_L*p02A + f_x*p12A - 2*f_x*p02A $$ $$ \frac{dP10A}{dt} = k_A*p00A*AraC_2 - k_{-A}*p10A - 2*k_R*p10A*LacI_4 $$ $$ \ + k_{-R}*p11A - f_x*p10A + f_x*p11A $$ $$ \frac{dP11A}{dt} = k_A*p01A*AraC_2 - k_{-A}*p11A + 2*k_R*p10A*LacI_4 $$ $$ \ - k_{-R}*p11A - k_R*p11A*LacI_4 + 2*k_{-R}*p12A - f_x*p11A $$ $$ \ - f_x*p11A + 2*f_x*p12A $$ $$ \frac{dP12A}{dt} = k_A*p02A*AraC_2 - k_{-A}*p12A + k_R*p11A*LacI_4 $$ $$ \ - 2*k_{-R}*p12A - k_L*p12A - f_x*p12A - 2*f_x*p12A $$ $$ \frac{dPL0A}{dt} = -k_UL*pL0A + 0.2*f_x*pL1A $$ $$ \frac{dPL1A}{dt} = 2*0.2*f_x*pL2A - 0.2*f_x*pL1A $$ $$ \frac{dPL2A}{dt} = k_L*p12A + k_L*p02A - 2*0.2*f_x*pL2A $$


Change in lacI promoter states over time:
$$ \frac{dP00L}{dt} = -k_A*p00L*AraC_2 + k_{-A}*p10L - 2*k_R*p00L*LacI_4 $$ $$ \ + k_{-R}*p01L + k_UL*pL0L + f_x*p10L + f_x*p01L $$ $$ \frac{dP01L}{dt} = -k_A*p01L*AraC_2 + k_{-A}*p11L + 2*k_R*p00L*LacI_4 $$ $$ \ - k_{-R}*p01L - k_R*p01L*LacI_4 + 2*k_{-R}*p02L + f_x*p11L $$ $$ \ - f_x*p01L + 2*f_x*p02L $$ $$ \frac{dP02L}{dt} = -k_A*p02L*AraC_2 + k_{-A}*p12L + k_R*p01L*LacI_4 $$ $$ \ - 2*k_{-R}*p02L - k_L*p02L + f_x*p12L - 2*f_x*p02L $$ $$ \frac{dP10L}{dt} = k_A*p00L*AraC_2 - k_{-A}*p10L - 2*k_R*p10L*LacI_4 $$ $$ \ + k_{-R}*p11L - f_x*p10L + f_x*p11L $$ $$ \frac{dP11L}{dt} = k_A*p01L*AraC_2 - k_{-A}*p11L + 2*k_R*p10L*LacI_4 $$ $$ \ - k_{-R}*p11L - k_R*p11L*LacI_4 + 2*k_{-R}*p12L - f_x*p11L $$ $$ \ - f_x*p11L + 2*f_x*p12L $$ $$ \frac{dP12L}{dt} = k_A*p02L*AraC_2 - k_{-A}*p12L + k_R*p11L*LacI_4 $$ $$ \ - 2*k_{-R}*p12L - k_L*p12L - f_x*p12L - 2*f_x*p12L $$ $$ \frac{dPL0L}{dt} = -k_UL*pL0L + 0.2*f_x*pL1L $$ $$ \frac{dPL1L}{dt} = 2*0.2*f_x*pL2L - 0.2*f_x*pL1L $$ $$ \frac{dGFP}{dt} = k_L*p12L + k_L*p02L - 2*0.2*f_x*pL2L $$


Change in gfp promoter states over time:
$$ \frac{dP00G}{dt} = -k_A*p00G*AraC_2 + k_{-A}*p10G - 2*k_R*p00G*LacI_4 $$ $$ \ + k_{-R}*p01G + k_UL*pL0G + f_x*p10G + f_x*p01G $$ $$ \frac{dP01G}{dt} = -k_A*p01G*AraC_2 + k_{-A}*p11G + 2*k_R*p00G*LacI_4 $$ $$ \ - k_{-R}*p01G - k_R*p01G*LacI_4 + 2*k_{-R}*p02G + f_x*p11G $$ $$ \ - f_x*p01G + 2*f_x*p02G $$ $$ \frac{dP02G}{dt} = -k_A*p02G*AraC_2 + k_{-A}*p12G + k_R*p01G*LacI_4 $$ $$ \ - 2*k_{-R}*p02G - k_L*p02G + f_x*p12G - 2*f_x*p02G $$ $$ \frac{dP10G}{dt} = k_A*p00G*AraC_2 - k_{-A}*p10G - 2*k_R*p10G*LacI_4 $$ $$ \ + k_{-R}*p11G - f_x*p10G + f_x*p11G $$ $$ \frac{dP11G}{dt} = k_A*p01G*AraC_2 - k_{-A}*p11G + 2*k_R*p10G*LacI_4 $$ $$ \ - k_{-R}*p11G - k_R*p11G*LacI_4 + 2*k_{-R}*p12G - f_x*p11G $$ $$ \ - f_x*p11G + 2*f_x*p12G $$ $$ \frac{dP12G}{dt} = k_A*p02G*AraC_2 - k_{-A}*p12G + k_R*p11G*LacI_4 $$ $$ \ - 2*k_{-R}*p12G - k_L*p12G - f_x*p12G - 2*f_x*p12GP $$ $$ \frac{dPL0G}{dt} = -k_UL*pL0G + 0.2*f_x*pL1G $$ $$ \frac{dPL1G}{dt} = 2*0.2*f_x*pL2G - 0.2*f_x*pL1G $$ $$ \frac{dPL2G}{dt} = k_L*p12G + k_L*p02G - 2*0.2*f_x*pL2G $$


Change in nif promoter states over time:
$$ \frac{dP00N}{dt} = -k_A*p00N*AraC_2 + k_{-A}*p10N - 2*k_R*p00N*LacI_4 $$ $$ \ + k_{-R}*p01N + k_UL*pL0N + f_x*p10N + f_x*p01N $$ $$ \frac{dP01N}{dt} = -k_A*p01N*AraC_2 + k_{-A}*p11N + 2*k_R*p00N*LacI_4 $$ $$ \ - k_{-R}*p01N - k_R*p01N*LacI_4 + 2*k_{-R}*p02N + f_x*p11N $$ $$ \ - f_x*p01N + 2*f_x*p02N $$ $$ \frac{dP02N}{dt} = -k_A*p02N*AraC_2 + k_{-A}*p12N + k_R*p01N*LacI_4 $$ $$ \ - 2*k_{-R}*p02N - k_L*p02N + f_x*p12N - 2*f_x*p02N $$ $$ \frac{dP10N}{dt} = k_A*p00N*AraC_2 - k_{-A}*p10N - 2*k_R*p10N*LacI_4 $$ $$ \ + k_{-R}*p11N - f_x*p10N + f_x*p11N $$ $$ \frac{dP11N}{dt} = k_A*p01N*AraC_2 - k_{-A}*p11N + 2*k_R*p10N*LacI_4 $$ $$ \ - k_{-R}*p11N - k_R*p11N*LacI_4 + 2*k_{-R}*p12N - f_x*p11N $$ $$ \ - f_x*p11N + 2*f_x*p12N $$ $$ \frac{dP12N}{dt} = k_A*p02N*AraC_2 - k_{-A}*p12N + k_R*p11N*LacI_4 $$ $$ \ - 2*k_{-R}*p12N - k_L*p12N - f_x*p12N - 2*f_x*p12N $$ $$ \frac{dPL0N}{dt} = -k_UL*pL0N + 0.2*f_x*pL1N $$ $$ \frac{dPL1N}{dt} = 2*0.2*f_x*pL2N - 0.2*f_x*pL1N $$ $$ \frac{dPL2N}{dt} = k_L*p12N + k_L*p02N - 2*0.2*f_x*pL2N $$






Nitrogenase Model

Scheme C is the flowchart of the Nitrogenase mechanism we modeled, incorporating negative cooperativity [3].


Assumptions

  • Fe protein does not bind to ADP

    • In reality, it does but it has lower affinity so it is negligible.

  • Re-reduction must occur before nucleotide exchange [5].

  • Rate of nucleotide exchange is the same as the rate of ATP binding to Fe [7].

  • Simultaneous binding of both Fe proteins to MoFe

    • Step 1: MoFe + Fe(ATP)2 --> MoFe(Fe(ATP)2), k (slow)

    • Step 2: MoFe(Fe(ATP)2) + Fe(ATP)2 --> MoFe(Fe(ATP)2)2 (fast)
  • Excess of [Fe4S4] clusters

    • This is not the limiting step.

  • Assumes excess FeMo-co and relatively fast insertion of FeMo-co compared to P-cluster maturation rate

  • ATP is assumed to be in excess so it's not a limiting factor

  • Excess NifX so transportation of M-cluster is not a limiting factor

  • Assembly of FeMo-co is assumed to be fast relative to p-cluster maturation

    • Excess FeMo-co compared to MoFe

  • MoFe is only a functional enzyme when 2 Fe proteins are bound to it

  • Fe protein formation is fast compared to NifH dimerization

  • Fe protein after formation is in its oxidized state

  • Fe protein must be reduced before it can bind to ATP

Nitrogenase code ran for 60 minutes with an initial concentration of 10000 molecules. Ammonia production was plotted on the y axis.


The above figure shows production of 20,000 molecules of NH3 in about 1 second. This matches literature that states a turnover rate of 1 N2 per second [4]. Since each molecule of N2 produces 2 NH3 molecules, 10,000 molecules of MoFe should produce 20,000 molecules of NH3 in 1 second. The subsequent decrease in slope accounts for the lower dissociation time using DT as a reactant versus an in vivo reductant like NifF, and the use of other various late limiting reactions [6]. Literature also suggests that the rates used in this model are accurate.


Nitrogenase Model— Nomenclature

FeOx Concentration of Fe protein in its oxidized state
FeRe Concentration of Fe protein in its reduced state
FeReATP2 Concentration of Fe protein in its reduced state and bound to 2 ATP molecules
FeOxADP2 Concentration of Fe protein in its oxidized state and bound to 2 ADP molecules
FeReADP2 Concentration of Fe protein in its reduced state and bound to 2 ADP molecules
FeReATP2MoFe0FeReATP2 Concentration of MoFe protein in its native state bound to 2 FeReATP2 complexes
FeReATP2MoFe1FeReATP2 Concentration of MoFe protein reduced with 1 electron bound to 2 FeReATP2 complexes
FeReATP2MoFe2FeReATP2 Concentration of MoFe protein reduced with 2 electron bound to 2 FeReATP2 complexes
FeReATP2MoFe3FeReATP2 Concentration of MoFe protein reduced with 3 electron bound to 2 FeReATP2 complexes
FeReATP2MoFe4FeReATP2 Concentration of MoFe protein reduced with 4 electron bound to 2 FeReATP2 complexes
FeReATP2MoFe5FeReATP2 Concentration of MoFe protein reduced with 5 electron bound to 2 FeReATP2 complexes
FeReATP2MoFe6FeReATP2 Concentration of MoFe protein reduced with 6 electron bound to 2 FeReATP2 complexes
FeReATP2MoFe7FeReATP2 MoFe protein reduced with 7 electron bound to 2 FeReATP2 complexes
MoFe0 Concentration of MoFe protein with no Fe proteins bound
MoFe2 Concentration of MoFe protein with no Fe proteins and 2 electron bound
MoFe3 Concentration of MoFe protein with no Fe proteins and 3 electron bound
MoFe4 Concentration of MoFe protein with no Fe proteins and 4 electron bound
MoFe5 Concentration of MoFe protein with no Fe proteins and 5 electron bound
MoFe6 Concentration of MoFe protein with no Fe proteins and 6 electron bound
MoFe7 Concentration of MoFe protein with no Fe proteins and 7 electron bound
MoFe0FeReATP2 Concentration of MoFe protein with 0 electrons and 1 FeReATP2 bound
MoFe1FeReATP2 Concentration of MoFe protein with 1 electrons and 1 FeReATP2 bound
MoFe2FeReATP2 Concentration of MoFe protein with 2 electrons and 1 FeReATP2 bound
MoFe3FeReATP2 Concentration of MoFe protein with 3 electrons and 1 FeReATP2 bound
MoFe4FeReATP2 Concentration of MoFe protein with 4 electrons and 1 FeReATP2 bound
MoFe5FeReATP2 Concentration of MoFe protein with 5 electrons and 1 FeReATP2 bound
MoFe6FeReATP2 Concentration of MoFe protein with 6 electrons and 1 FeReATP2 bound
MoFe7FeReATP2 Concentration of MoFe protein with 7 electrons and 1 FeReATP2 bound
FeOxATP2MoFe1FeReATP2 Concentration of MoFe protein with 1 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe2FeReATP2 Concentration of MoFe protein with 2 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe3FeReATP2 Concentration of MoFe protein with 3 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe4FeReATP2 Concentration of MoFe protein with 4 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe5FeReATP2 Concentration of MoFe protein with 5 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe6FeReATP2 Concentration of MoFe protein with 6 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe7FeReATP2 Concentration of MoFe protein with 7 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxATP2MoFe8FeReATP2 Concentration of MoFe protein with 8 electrons and 2 Fe proteins bound after 1 electron transfer
FeOxADPPi2MoFe1FeReATP2 Concentration of MoFe protein with 1 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe2FeReATP2 Concentration of MoFe protein with 2 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe3FeReATP2 Concentration of MoFe protein with 3 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe4FeReATP2 Concentration of MoFe protein with 4 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe5FeReATP2 Concentration of MoFe protein with 5 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe6FeReATP2 Concentration of MoFe protein with 6 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe7FeReATP2 Concentration of MoFe protein with 7 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe8FeReATP2 Concentration of MoFe protein with 8 electrons and 2 Fe proteins bound after ATP hydrolysis
FeOxADP2MoFe1FeReATP2 Concentration of MoFe protein with 1 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe2FeReATP2 Concentration of MoFe protein with 2 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe3FeReATP2 Concentration of MoFe protein with 3 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe4FeReATP2 Concentration of MoFe protein with 4 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe5FeReATP2 Concentration of MoFe protein with 5 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe6FeReATP2 Concentration of MoFe protein with 6 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe7FeReATP2 Concentration of MoFe protein with 7 electrons and 2 Fe proteins bound after phosphate release
FeOxADP2MoFe8FeReATP2 Concentration of MoFe protein with 8 electrons and 2 Fe proteins bound after phosphate release
FeOxATP2MoFe2FeOxATP2 Concentration of MoFe protein with 2 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxATP2MoFe3FeOxATP2 Concentration of MoFe protein with 3 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxATP2MoFe4FeOxATP2 Concentration of MoFe protein with 4 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxATP2MoFe5FeOxATP2 Concentration of MoFe protein with 5 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxATP2MoFe6FeOxATP2 Concentration of MoFe protein with 6 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxATP2MoFe7FeOxATP2 Concentration of MoFe protein with 7 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxATP2MoFe8FeOxATP2 Concentration of MoFe protein with 8 electrons and 2 Fe proteins bound after transfer of electrons from both Fe proteins
FeOxADPPi2MoFe2FeOxADPPi2 Concentration of MoFe protein with 2 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe3FeOxADPPi2 Concentration of MoFe protein with 3 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe4FeOxADPPi2 Concentration of MoFe protein with 4 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe5FeOxADPPi2 Concentration of MoFe protein with 5 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe6FeOxADPPi2 Concentration of MoFe protein with 6 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe7FeOxADPPi2 Concentration of MoFe protein with 7 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADPPi2MoFe8FeOxADPPi2 Concentration of MoFe protein with 8 electrons and 2 oxidized Fe proteins bound after ATP hydrolysis
FeOxADP2MoFe2FeOxADP2 Concentration of MoFe protein with 2 electrons and 2 oxidized Fe proteins bound after phosphate release
FeOxADP2MoFe3FeOxADP2 Concentration of MoFe protein with 3 electrons and 2 oxidized Fe proteins bound after phosphate release
FeOxADP2MoFe4FeOxADP2 Concentration of MoFe protein with 4 electrons and 2 oxidized Fe proteins bound after phosphate release
FeOxADP2MoFe5FeOxADP2 Concentration of MoFe protein with 5 electrons and 2 oxidized Fe proteins bound after phosphate release
FeOxADP2MoFe6FeOxADP2 Concentration of MoFe protein with 6 electrons and 2 oxidized Fe proteins bound after phosphate release
FeOxADP2MoFe7FeOxADP2 Concentration of MoFe protein with 7 electrons and 2 oxidized Fe proteins bound after phosphate release
FeOxADP2MoFe8FeOxADP2 Concentration of MoFe protein with 8 electrons and 2 oxidized Fe proteins bound after phosphate release
NH3 Ammonia Concentration

Nitrogenase Model— Parameters

Parameter (s-1) Description
kFb 0.033211 [41] Rate of binding of Fere(ATP)2 to MoFe protein
kEt 267 ± 20 [43] Rate of electron transfer from 1st Fe protein to MoFe protein
kH1 62 ± 4 [43] Rate of ATP hydrolysis from 1st Fe protein
kPr1 19 ± 1 [43] Rate of phosphate release from 1st Fe protein
kD1 13 ± 2 [43] Rate of Feox(ADP)2 dissociation from MoFe complex
kEt2 37 ± 6 [43] Rate of electron transfer from 2nd Fe protein
kH2 31 ± 3 [43] Rate of ATP hydrolysis from 2nd Fe protein
kPr2 0.6 ± 5 [43] Rate of phosphate release from 2nd Fe protein
kD2 30 ± 25 [43] Rate of both Feox(ADP)2 dissociation together from MoFe
kFare 1538 [42] Rate of Feox(ADP)2 reduction
kNs 70 [44] Rate of nucleotide swap from Fere(ADP)2 to Feox(ATP)2
kFre 1188 [42] Rate of unbound Fe protein reduction

MATLAB code for our Nitrogenase Model equations is embedded below.


dMoFe3FeReATP2 = kD1*FeOxADP2MoFe3FeReATP2 -kFb*MoFe3FeReATP2*FeReATP2 + kFb*MoFe3*FeReATP2 dMoFe4FeReATP2 = kD1*FeOxADP2MoFe4FeReATP2 -kFb*MoFe4FeReATP2*FeReATP2 + kFb*MoFe4*FeReATP2 dMoFe5FeReATP2 = kD1*FeOxADP2MoFe5FeReATP2 -kFb*MoFe5FeReATP2*FeReATP2 + kFb*MoFe5*FeReATP2 dMoFe6FeReATP2 = kD1*FeOxADP2MoFe6FeReATP2 -kFb*MoFe6FeReATP2*FeReATP2 + kFb*MoFe6*FeReATP2 dMoFe7FeReATP2 = kD1*FeOxADP2MoFe7FeReATP2 -kFb*MoFe7FeReATP2*FeReATP2 + kFb*MoFe7*FeReATP2


dFeOxATP2MoFe1FeReATP2 = kEt*FeReATP2MoFe0FeReATP2 - kH1*FeOxATP2MoFe1FeReATP2 - kEt2*FeOxATP2MoFe1FeReATP2 dFeOxATP2MoFe2FeReATP2 = kEt*FeReATP2MoFe1FeReATP2 - kH1*FeOxATP2MoFe2FeReATP2 - kEt2*FeOxATP2MoFe2FeReATP2 dFeOxATP2MoFe3FeReATP2 = kEt*FeReATP2MoFe2FeReATP2 - kH1*FeOxATP2MoFe3FeReATP2 - kEt2*FeOxATP2MoFe3FeReATP2 dFeOxATP2MoFe4FeReATP2 = kEt*FeReATP2MoFe3FeReATP2 - kH1*FeOxATP2MoFe4FeReATP2 - kEt2*FeOxATP2MoFe4FeReATP2 dFeOxATP2MoFe5FeReATP2 = kEt*FeReATP2MoFe4FeReATP2 - kH1*FeOxATP2MoFe5FeReATP2 - kEt2*FeOxATP2MoFe5FeReATP2 dFeOxATP2MoFe6FeReATP2 = kEt*FeReATP2MoFe5FeReATP2 - kH1*FeOxATP2MoFe6FeReATP2 - kEt2*FeOxATP2MoFe6FeReATP2 dFeOxATP2MoFe7FeReATP2 = kEt*FeReATP2MoFe6FeReATP2 - kH1*FeOxATP2MoFe7FeReATP2 - kEt2*FeOxATP2MoFe7FeReATP2 dFeOxATP2MoFe8FeReATP2 = kEt*FeReATP2MoFe7FeReATP2 - kH1*FeOxATP2MoFe8FeReATP2


dFeOxADPPi2MoFe1FeReATP2 = kH1*FeOxATP2MoFe1FeReATP2 - kPr1*FeOxADPPi2MoFe1FeReATP2 dFeOxADPPi2MoFe2FeReATP2 = kH1*FeOxATP2MoFe2FeReATP2 - kPr1*FeOxADPPi2MoFe2FeReATP2 dFeOxADPPi2MoFe3FeReATP2 = kH1*FeOxATP2MoFe3FeReATP2 - kPr1*FeOxADPPi2MoFe3FeReATP2 dFeOxADPPi2MoFe4FeReATP2 = kH1*FeOxATP2MoFe4FeReATP2 - kPr1*FeOxADPPi2MoFe4FeReATP2 dFeOxADPPi2MoFe5FeReATP2 = kH1*FeOxATP2MoFe5FeReATP2 - kPr1*FeOxADPPi2MoFe5FeReATP2 dFeOxADPPi2MoFe6FeReATP2 = kH1*FeOxATP2MoFe6FeReATP2 - kPr1*FeOxADPPi2MoFe6FeReATP2 dFeOxADPPi2MoFe7FeReATP2 = kH1*FeOxATP2MoFe7FeReATP2 - kPr1*FeOxADPPi2MoFe7FeReATP2 dFeOxADPPi2MoFe8FeReATP2 = kH1*FeOxATP2MoFe8FeReATP2 - kPr1*FeOxADPPi2MoFe8FeReATP2


dFeOxADP2MoFe1FeReATP2 = kPr1*FeOxADPPi2MoFe1FeReATP2 - kD1*FeOxADP2MoFe1FeReATP2 dFeOxADP2MoFe2FeReATP2 = kPr1*FeOxADPPi2MoFe2FeReATP2 - kD1*FeOxADP2MoFe2FeReATP2 dFeOxADP2MoFe3FeReATP2 = kPr1*FeOxADPPi2MoFe3FeReATP2 - kD1*FeOxADP2MoFe3FeReATP2 dFeOxADP2MoFe4FeReATP2 = kPr1*FeOxADPPi2MoFe4FeReATP2 - kD1*FeOxADP2MoFe4FeReATP2 dFeOxADP2MoFe5FeReATP2 = kPr1*FeOxADPPi2MoFe5FeReATP2 - kD1*FeOxADP2MoFe5FeReATP2 dFeOxADP2MoFe6FeReATP2 = kPr1*FeOxADPPi2MoFe6FeReATP2 - kD1*FeOxADP2MoFe6FeReATP2 dFeOxADP2MoFe7FeReATP2 = kPr1*FeOxADPPi2MoFe7FeReATP2 - kD1*FeOxADP2MoFe7FeReATP2 dFeOxADP2MoFe8FeReATP2 = kPr1*FeOxADPPi2MoFe8FeReATP2 - kD1*FeOxADP2MoFe8FeReATP2


dFeOxATP2MoFe2FeOxATP2 = kEt2*FeOxATP2MoFe1FeReATP2 - kH2*FeOxATP2MoFe2FeOxATP2 dFeOxATP2MoFe3FeOxATP2 = kEt2*FeOxATP2MoFe2FeReATP2 - kH2*FeOxATP2MoFe3FeOxATP2 dFeOxATP2MoFe4FeOxATP2 = kEt2*FeOxATP2MoFe3FeReATP2 - kH2*FeOxATP2MoFe4FeOxATP2 dFeOxATP2MoFe5FeOxATP2 = kEt2*FeOxATP2MoFe4FeReATP2 - kH2*FeOxATP2MoFe5FeOxATP2 dFeOxATP2MoFe6FeOxATP2 = kEt2*FeOxATP2MoFe5FeReATP2 - kH2*FeOxATP2MoFe6FeOxATP2 dFeOxATP2MoFe7FeOxATP2 = kEt2*FeOxATP2MoFe6FeReATP2 - kH2*FeOxATP2MoFe7FeOxATP2 dFeOxATP2MoFe8FeOxATP2 = kEt2*FeOxATP2MoFe7FeReATP2 - kH2*FeOxATP2MoFe8FeOxATP2


dFeOxADPPi2MoFe2FeOxADPPi2 = kH2*FeOxATP2MoFe2FeOxATP2 - kPr2*FeOxADPPi2MoFe2FeOxADPPi2 dFeOxADPPi2MoFe3FeOxADPPi2 = kH2*FeOxATP2MoFe3FeOxATP2 - kPr2*FeOxADPPi2MoFe3FeOxADPPi2 dFeOxADPPi2MoFe4FeOxADPPi2 = kH2*FeOxATP2MoFe4FeOxATP2 - kPr2*FeOxADPPi2MoFe4FeOxADPPi2 dFeOxADPPi2MoFe5FeOxADPPi2 = kH2*FeOxATP2MoFe5FeOxATP2 - kPr2*FeOxADPPi2MoFe5FeOxADPPi2 dFeOxADPPi2MoFe6FeOxADPPi2 = kH2*FeOxATP2MoFe6FeOxATP2 - kPr2*FeOxADPPi2MoFe6FeOxADPPi2 dFeOxADPPi2MoFe7FeOxADPPi2 = kH2*FeOxATP2MoFe7FeOxATP2 - kPr2*FeOxADPPi2MoFe7FeOxADPPi2 dFeOxADPPi2MoFe8FeOxADPPi2 = kH2*FeOxATP2MoFe8FeOxATP2 - kPr2*FeOxADPPi2MoFe8FeOxADPPi2


dFeOxADP2MoFe2FeOxADP2 = kPr2*FeOxADPPi2MoFe2FeOxADPPi2 - kD2*FeOxADP2MoFe2FeOxADP2 dFeOxADP2MoFe3FeOxADP2 = kPr2*FeOxADPPi2MoFe3FeOxADPPi2 - kD2*FeOxADP2MoFe3FeOxADP2 dFeOxADP2MoFe4FeOxADP2 = kPr2*FeOxADPPi2MoFe4FeOxADPPi2 - kD2*FeOxADP2MoFe4FeOxADP2 dFeOxADP2MoFe5FeOxADP2 = kPr2*FeOxADPPi2MoFe5FeOxADPPi2 - kD2*FeOxADP2MoFe5FeOxADP2 dFeOxADP2MoFe6FeOxADP2 = kPr2*FeOxADPPi2MoFe6FeOxADPPi2 - kD2*FeOxADP2MoFe6FeOxADP2 dFeOxADP2MoFe7FeOxADP2 = kPr2*FeOxADPPi2MoFe7FeOxADPPi2 - kD2*FeOxADP2MoFe7FeOxADP2 dFeOxADP2MoFe8FeOxADP2 = kPr2*FeOxADPPi2MoFe8FeOxADPPi2 - kD2*FeOxADP2MoFe8FeOxADP2


dNH3 = 2 * (kD1*FeOxADP2MoFe8FeReATP2 + 2*kD2*FeOxADP2MoFe8FeOxADP2)


Nitrogenase code that cuts off the ability to regenerate Fe proteins after the reduction of N2.


To verify the lack of errors, the code cut off the ability of nitrogenase to cycle at different periods of time. The ability of the code to have a steady state concentration at 20,000 molecules verifies the lack of bugs in the code. This was repeated for each of the subsequent steps of the mechanism, and we were able to consistently obtain desired concentrations of product, depending on which steady state we were observing.


Concentrations of ammonia run over a period of 32 minutes. The Nitroscillator code was run over a period of time that the genetic oscillator verified to contain 2 oscillations.


After verifying to the best of our ability that our nitrogenase catalysis rate matched literature and was accurate, we attempted to incorporate the nitrogenase code into the genetic oscillator. Our goal was to obtain oscillating changes in fixation rate, from periods of high fixation rates to periods of little to no fixation rates. This was verified by observing ammonia concentrations. Due to the complexity of the code, we chose the arabinose and IPTG conditions that gave the shortest periods deducted from the table we created. The above figure shows two oscillations in nitrogen fixation rates, demonstrated by a spike in protein expression followed by a flat steady state period.





References

  1. Jiang, X., Payá-Tormo, L., Coroian, D. et al. Exploiting genetic diversity and gene synthesis to identify superior nitrogenase NifH protein variants to engineer N2-fixation in plants. Commun Biol 4, 4 (2021). https://doi.org/10.1038/s42003-020-01536-6

  2. Stricker, J., Cookson, S., Bennett, M. R., Mather, W. H., Tsimring, L. S., & Hasty, J. (2008). A fast, robust and tunable synthetic gene oscillator. Nature, 456(7221), 516–519. https://doi.org/10.1038/nature07389

  3. Danyal, Karamatullah, et al. “Negative Cooperativity in the Nitrogenase Fe Protein Electron Delivery Cycle.” Proceedings of the National Academy of Sciences of the United States of America, vol. 113, no. 40, 3 Oct. 2016, https://doi.org/10.1073/pnas.1613089113.

  4. Einsle, Oliver, and Douglas C. Rees. “Structural Enzymology of Nitrogenase Enzymes.” Chemical Reviews, vol. 120, no. 12, 15 June 2020, pp. 4969–5004, https://doi.org/10.1021/acs.chemrev.0c00067.

  5. Pence, Natasha, et al. “Unraveling the Interactions of the Physiological Reductant Flavodoxin with the Different Conformations of the Fe Protein in the Nitrogenase Cycle.” Biological Chemistry, vol. 292, no. 38, 22 Sept. 2017, pp. 15661–15669, https://doi.org/10.1074/jbc.m117.801548.

  6. Rutledge, Hannah L., and F. Akif Tezcan. “Electron Transfer in Nitrogenase.” Chemical Reviews, vol. 120, no. 12, 30 Jan. 2020, pp. 5158–5193, www.ncbi.nlm.nih.gov/pmc/articles/PMC7466952/, https://doi.org/10.1021/acs.chemrev.9b00663.

  7. Shah, Sudipta. Role of ATP Hydrolysis and Mechanism of Substrate Reduction in Nitrogenase. May 2017.