Part I: Fitting predicted
Bacillus subtilis
rowth curves and exploring the optimal conditions for its growth environment
(1)Exploring the effect of temperature on the growth of
Bacillus subtilis
(pH=7). OD
630
on the vertical Y-axis against time (hours) on the horizontal X-axis.
Figure 1.
Bacillus subtilis
rowth curve at 16°C (100ml of medium inoculated with 0.4ml of bacterial solution)
Figure 2.
Bacillus subtilis
growth curve at 25°C (100ml of medium inoculated with 0.4ml of bacterial solution)
At lower temperatures, the growth curve of
Bacillus subtilis
experienced a flattening out followed by a sharp rise, showing that
Bacillus subtilis
grows more slowly at certain lower temperatures.
Figure 3.
Bacillus subtilis
growth curve at 36°C (100ml of medium inoculated with 0.4ml of bacterial solution)
Figure 4.
Finding the K value and fitting the equation
The growth curve of
Bacillus subtilis
at this temperature showed a smooth and then increasing and then decreasing curve change, and it is expected that the OD value measured at this temperature reached 2.068.
The analysis of the three graphs shows that the higher the temperature the faster the growth rate of
Bacillus subtilis
within a certain temperature range. It is evident that the higher the temperature the more suitable the growth of
Bacillus subtilis
within a certain temperature range.
(2)Exploring the effect of pH on the growth of
Bacillus subtilis
(at 36°C)
Figure 5.
Effect of different initial pH values on the growth of the strain
Under acidic or alkaline conditions,
Bacillus subtilis
may produce budding spores that are not suitable for the experimental conditions. After consulting relevant papers and materials, we concluded that at an initial pH of 6.0-7.2, the OD
600
value of
Bacillus subtilis
fermentation broth increases with the initial pH, and the OD
600
value of the strain tends to decrease when the pH is greater than 7.2. Therefore, the optimum initial pH of the
Bacillus subtilis
strain is 7.2.
Experimental data as a guide to real life
1. When we use
Bacillus subtilis
in fish tanks, we should monitor the temperature of the water in the tank through the temperature sensor in real time within a reasonable temperature range, so that the temperature is close to the optimum growth temperature of
Bacillus subtilis
can achieve the fastest growth rate and thus achieve the best effect of tank treatment.
2.If the pH is too low and the water is acidic, put some coral sand in the tank or add the right amount of sodium bicarbonate; if the pH is too high and the water is alkaline, add olive leaves, pineapple bark, sunken wood, etc. or add the right amount of dihydrogen carbonate to lower the pH and make the water neutral, so that the
Bacillus subtilis
can rise quickly. can be quickly upgraded to achieve the effect of rapid treatment of sewage.
Part II: Using SimBiology plates to explore the quantitative changes of substances in the metabolic pathway of heterotrophic nitrification- aerobic denitrification denitrification
(i) Determining chemical reaction equations
Figure 6.
Heterotrophic nitrification - aerobic denitrification metabolic pathway for denitrification
Figure 7.
Principle of heterotrophic nitrification - aerobic denitrification denitrification metabolism
Based on the heterotrophic nitrification - aerobic denitrification denitrification reaction, we can derive the following chemical reaction equation by consulting the literature, etc.
3 NH
3
+ 0
2
<-> 2NH
2
OH + 2 H
2
0; (2.1)
2 NH
3
OH + 2 H
2
0 <-> 2 NO
2
-
; (2.2)
2 NO
2
-
<-> 2 NO; (2.3)
2 NO <-> N
2
0; (2.4)
N
2
0 <-> N
2; (2.5)
(ii) Modelling with SimBiology
Figure 8.
Flow of the biological reaction process modelled and simulated by SimBiology
Based on this list of five chemical reaction equations, we determine the relationship between the masses of the substances and the reaction equation, so that the changes in the substances are linked to each other, according to this reaction equation, we use the simbiology module of MATLAB to model it first.
1. Open the SimBi-ology interface by typing "sbiodesktop" in the MATLAB command window. In the "File/NewProject/ProjectType" list box, select "Model" and click "OK". Double click on the "ProjectExplorer" on the left side of the interface and click on the "File/SaveProjectas" option and type in the name of the project file and click "Save ".
2. Double click on 'ProjectExplorer/TableView/Reactions' and type the reaction equation in the 'Enter Reaction' box.
3. Build the model in Simbiology using "Rules": Enter the differential rate equation or algebraic equation in the "Rules" editor window.
4. Modeling with "KineticLaw" in Simbiology: Select the kinetic model in the "KineticLaw" list or customize it.
5. Click on "New" in the "SpecifyCorespondingParameterNames" box and type in the "Name In the "Name" box type in the name of the parameter and in the "Scope" column select "KineticLaw" and click on "Create".
6. In the "ModelSesion-untitle" dialog box click on "Species" in the "InitialAmount In the "InitialAmount" cell box and the "InitialAmountUnits" list, type in the initial values for the products and reactants and select the units for them to complete the initial value setting for the substance.
7. At any time during the modelling process, click on "Ver-ify" at the top right of the SimBiology screen to verify that the model and parameters entered are not incorrect.
8. Use the menu "Analysis/AddAnalysisTasktountitled/Simu-latemodel" to open the "Simulation" box and set the simulation conditions. In the upper right-hand corner, click on "R un (Run)" to run the simulation of the input model.
Figure 9.
Modelling of the reaction model in simbiology
By creating five plates A, B, C, D and E so that they correspond to five reaction equations, by creating NH
2,
O
2,
NH
2
OH, H
2
O, NO
2 -,
NO, N
2
O
,N
2
substances within each of the five plates A,B,C,D,E, followed by determining their initial amount of substance, assuming NH
3
is 3 mol, 0
2
is 1 mol, H
2
0 is 1 mol and the rest of the The initial amount of substance for each of the remaining substances is 0 mol.
After determining the amount of substances and the connection between the individual reaction equations, we used simbiology simulation runs to analyse how each substance changed during this reaction.
Figure 10.
Changes in the amount of material in denitrification by heterotrophic nitrification - aerobic denitrification
The graph of the change in the amount of substances shows the change in the content of NH
2
OH, NO and N
2
O toxic substances to a certain extent, so that measures can be taken in advance to prevent the impact of such toxic substances on the laboratory and laboratory personnel. The modelling and the amount of substances can be modified in real time according to the requirements of the experiment, thus helping the experiment to proceed better. By adjusting the amount of reactants, we can adjust the process of the reaction and detect whether the EL222-activated genes are properly expressed in real time by combining the target products AMO, HAO, NXR, NAR, NIR, NOR and NOS.
Introduction
In the initial stage of the experiment, we learned the common ways of wastewater nitrogen removal through literature reading, and determined that the nitrogen removal pathway adopted in this project was to express the following four proteases in
Bacillus subtilis,
ammonia monooxygenase (AMO), hydroxylamine oxidase (HAO), assimilatory nitroreductase (NAR), and nitroreductase (NIR)
to treat sewage. But the Research Collaboratory for Structural Bioinformatics (RCSB) doesn't have structural files for these proteins. Based on this situation, we used
Alphafold2
to predict the structures of four enzyme molecules. Subsequently, we used
molecular dynamics simulation (MDS)
to optimize the predicted structure of Alphafold2 to make it closer to the natural conformation and evaluated the stability of each region of the enzyme molecule to analyze the unstable regions in the enzyme molecule. In addition, we
predicted the region of the enzyme activity center
of each enzyme molecule through molecular docking, providing reference information for the structural design and optimization of four enzyme molecules.
Structure prediction
The amino acid sequences of the four enzyme molecules were obtained from the NCBI database
(Table 1)
, and Alphafold2 was used to predict the three-dimensional structure of the proteins from the amino acid sequences since the structure files of the corresponding enzyme molecules were not currently available.
Based on the pLDDT values, the 3D structures obtained from Alphafold2 were all considered to be of high quality, with the pLDDT value of HAO even reaching
98.24
(Table 1)
.
Table 1.
Information on the four proteins in the experiment.
The following Figures 1 to 4 show the 3D models of AMO, HAO, NAR, and NIR obtained using Alphafold2, respectively.
Figure 1.
Predicted structure of AMO (pLDDT = 91.02).
Figure 2.
Predicted structure of HAO (pLDDT = 98.24).
Figure 3.
Predicted structure of NAR (pLDDT = 88.09).
Figure 4.
Predicted structure of NIR (pLDDT = 92.83).
Docking
We
predicted the pockets
of the enzyme molecules after acquiring the appropriate 3D structure data, and the prediction results show that HAO and NIR molecules correctly predicted the high-quality pockets
(Figure 5, 6)
, but AMO and NAR did not acquire suitable pockets
(Table 3)
.
As a result,
we performed molecular docking for HAO and NIR, resulting in high-quality pockets with their matching substrates.
We intend to repeat pocket prediction for the resultant AMO and NAR conformations following molecular dynamics simulations for AMO and NAR because the pocket of some proteins may be revealed after conformational optimization using molecular dynamics simulations.
Table 2.
Pocket scores and probability of the structures predicted by Alphafold2.
Figure 5.
3D display of the predicted HAO pocket (red).
Figure 6.
3D display of the predicted NIR pocket (blue).
Based on the molecular docking data, we investigated the interactions between the substrate and the enzyme molecules, and the results are presented in the pictures below.
The molecular docking of HAO with NH
2
OH revealed that the lowest binding energy of the NH
2
OH molecule was
-2.8 kcal/mol
and four amino acid residues
(His 188, Glu 190, Asn 198, Arg 201)
in the active center of HAO were able to form
six hydrogen bonds
with the best conformation of the NH
2
OH molecule; Furthermore, three amino acid residues
(Met 189, Tyr 191, and Ile 199)
developed hydrophobic interactions with the hydroxylamine molecule
(Figure7, 8).
Figure 7.
HAO amino acid residues react with NH
2
OH in three dimensions.
Figure 8.
HAO amino acid residues react with NH
2
OH in two dimensions.
The docking of NIR with the NO
2 -
on reveals that the NO
2-
makes a hydrogen bond with each of the three amino acid residues
(Arg 399, Lys 454, and Glu 512)
in the NIR pocket, and they all interact with the NO
2-
oxygen atom. The tyrosine at position 570 forms
hydrophobic interactions
with the nicotinic acid ion. Interaction
(Figure 9, 10)
The best binding affinity to the nitrite ion through NIR is -3.4 kcal/mol.
Figure 9.
NIR amino acid residues react with NO
2 -
in two dimensions.
Figure 10.
NIR amino acid residues react with NO
2 -
in two dimensions.
Molecular dynamics simulation
We next ran
100 ns molecular dynamics simulations
(MDS) to improve the structure of the acquired enzyme molecule and get it closer to its native shape. Simultaneously, we hope to evaluate the stability of each part of the enzyme molecule using molecular dynamics simulation to identify the flexible regions in the enzyme molecule and provide a theoretical reference for subsequent modification and optimization of the enzyme molecule to
improve its stability and thus its catalytic efficiency.
Furthermore, as previously stated, molecular dynamics simulations may
bring out the hidden active centers
in some enzyme molecules, and we intend to further investigate the position of the pockets of AMO and NAR in the simulated conformations to study their interactions with substrate molecules.
We detailed the
GROMACS
analysis procedure of GXU-China in 2022 by referring to the GROMACS tutorial and making relevant parameter modifications. In summary, we utilized GROMACS to construct the protease structure's topology file, then placed it in a box, filled it with water solvent, and added ions to balance the protein's charge. The system was then subjected to
energy minimization,
NVT equilibrium,
and
NPT equilibrium,
after which we performed 100 ns MDS and presented the protein's movement track in the simulation process.
The MDS results of the four enzyme molecules chosen for the study are as follows:
First, we visualized the trajectories of the four enzyme molecules' molecular dynamics simulations
(Figure 11)
, in order to view and assess the conformational changes of the enzyme molecules during the simulations.
(A)
AMO
(B)
HAO
(C)
NAR
(D)
NIR
Figure 11.
Animation of the molecular dynamics simulation trajectories of the four enzyme molecules (part).
After 100 ns of simulation, the conformations of the four enzyme molecules were aligned to the pre-simulation conformation, and the
Root Mean Square Deviation (RMSD)
was obtained. The pre-simulation conformation is shown in red, whereas the post-simulation conformation is highlighted in blue. The 100 ns of simulation makes the enzyme molecule conception more compatible with the natural conformation. According to the RMSD values,
AMO and NIR have less conformational changes and more stable structures before and after the simulation, but the other two enzyme molecules, NAR and HAO, have undergone significant conformational changes.
(A)
AMO (RMSD = 3.658 Å)
(B)
HAO (RMSD = 5.438 Å)
(C)
NAR (RMSD = 5.769 Å)
(D)
NIR (RMSD = 2.225 Å)
Figure 12.
Comparison of protein conformation before (red) and after (blue) 100 ns molecular dynamics simulation.
The simulation results file was used to compute and
show the changes in RMSD of the backbone of the four enzyme molecules during the experiment
(Figure 13)
. After 75 ns, the states of all four enzyme molecules essentially achieved equilibrium, and
AMO demonstrated the best stability
during the simulation with a
steady RMSD of about 4 Å
, followed by NIR, which was consistent with the results of the previous simulation's structural comparison.
Figure 13.
The RMSD of AMO, HAO, NAR, and NIR in backbone during 100 ns simulation.
We computed the
Root Mean Square Fluctuation (RMSF)
and
b-factor
of each structure from the result files to better understand the stability of the structures of each section of the four enzyme molecules and visualized the results as follows
(Figure 14-21)
.
We can observe from the distribution of RMSF and b-factor in each point of the AMO structure that
the AMO structure is fairly stable in the segment of about 1-125, with RMSF around 2 Å
. The remaining region, the region about 125–175, where the maximum RMSF exceeds 6 Å, is the least stable and contributes significantly to the overall instability of the AMO structure.
The model of b-factor coloring also shows that
this area is a significant portion of a loop and another part of a pair of β-fold sheets,
and future enzyme molecule modifications can target these regions by adding
proline(Pro) mutations
to improve the overall stability of the enzyme molecule
(Figure 14, 15)
.
Figure 14.
RMSF of each amino acid residue of AMO molecule during 100 ns molecular dynamics simulation.
Figure 15.
Structure of AMO colored by b-factor.
Similarly, the RMSF and b-factor distributions of HAO show that there are
three unstable areas
in the structure, centered around
0–50
, about
125
, and after
250
, with RMSF values of
250
, with RMSF values of
more than 5 Å
. The structural model colored according to the b-factor shows that these sections are the
uneven curl
of the protein molecule's C-and N-terminal, and
the loop structure in the middle position is therefore less stable
(Figure 16, 17)
.
Figure 16.
RMSF of each amino acid residue of HAO molecule during 100 ns molecular dynamics simulation.
Figure 17.
Structure of HAO colored by b-factor.
NAR's structure is
mainly stable
, with unstable pieces in smaller sections around 0, 200, 500, and 600, which are all unevenly coiled, as shown in the protein structure
(Figure 18, 19)
.
Figure 18.
RMSF of each amino acid residue of NAR molecule during 100 ns molecular dynamics simulation.
Figure 19.
Structure of NAR colored by b-factor.
In contrast to the three previously known protein structures,
NIR exhibits remarkable stability
, with the exception of the severely unstable region between 0-50. When combined with the b-factor coloring model, it is clear that the most unstable region is the uneven curl at the structure's end, but the other regions contain
huge segments of β-folded sheets
and
hence have great stability
(Figure 20, 21)
.
Figure 20.
RMSF of each amino acid residue of NIR molecule during 100 ns molecular dynamics simulation.
Figure 21.
Structure of NIR colored by b-factor.
Post-MDS docking
Previously, we discovered that the architectures of AMO and NAR derived using Alphafold2 could not be anticipated to produce high-quality pockets. Because the active centers of some enzyme molecules may be revealed after molecular dynamics simulation optimization, we
repeated pocket prediction using the protein conformation acquired
after molecular dynamics simulation.
Table 3.
Pocket scores and probability levels of AMO and NAR after molecular dynamics
Following the simulation, the projected results of the pockets of AMO and NAR conformations show that the scores and
probability of the predicted pockets have improved significantly
, and it is now feasible to
get high-quality pockets
that fulfill the requirements
(Figure 22, 23, Table 3)
.
This also supports our earlier hypothesis that
the active centers of AMO and NAR were discovered after 100 ns of conformational optimization using molecular dynamics modeling
.
Figure 22.
3D display of the predicted AMO pocket (green).
Figure 23.
3D display of the predicted NAR pocket (orange).
We performed molecular docking of AMO and NAR with their respective substrates using high-quality pockets and optimized conformations derived from molecular dynamics simulations and evaluated the enzyme-substrate interactions based on the docking results.
The NH
4+
establishes
hydrogen bonding contacts
with each of the four amino acid residues in the AMO pocket
(Trp 9, Trp 14, Trp 15, Thr 19)
and
hydrophobic interactions
with the
Val 18
. Their lowest binding affinity is
-1.4 kcal/mol
(Figure 24)
.
(A)
AMO amino acid residues react with NH
4+
in three dimensions.
(B)
AMO amino acid residues react with NH
4+
in two dimensions.
Figure 24.
Docking results of AMO and NH
4 +.
The NO
3 -
interacts with the amino acids in the NAR pocket exclusively through the
oxygen atom
, with the three oxygen atoms establishing
hydrogen bonds
with
Lys 541, Tyr 546, and Lys 583,
respectively, as well as
hydrophobic interactions
with
His 541 and Leu 586.
The binding energy was
-3.8 kcal/mol
,
the lowest of the four groups
of enzyme molecules and substrates studied
(Figure 25)
.
(A)
NAR amino acid residues react with NO
3 -
in three dimensions.
(B)
NAR amino acid residues react with NO
3 -
in two dimensions.
Figure 25.
Docking results of NAR and NO
3 -
.
Furthermore, we investigated the
binding statuses
of NIR and HAO to their respective substrates before and after the simulations. The scores and probability of NIR were both reduced significantly but remained reliable; the score values of HAO also declined slightly, but the confidence level did not change and remained at a high level
(Table 4).
Table 4.
Pocket scores and probabilities of NIR and HAO after MDS.
Similarly, we performed molecular docking with the substrate using the new conformation and the expected pocket information, examined the status of the substrate-enzyme molecule interaction at this moment, and compared the docking findings to those of the pre-simulation conformation.
The comparative results show that
the binding locations of substrate molecules in the pocket altered before and after the simulation in both NIR and HAO
(Figure 26 A and B).
The substrate molecule's binding energy to the enzyme molecule changed somewhat as well.
The binding affinity of the HAO molecule at the new optimum binding location rose by 0.2 kcal/mol to -2.6 kcal/mol, whereas the binding energy of NIR was reduced to -3.6 kcal/mol, improving stability marginally over the pre-simulation binding state.
The binding states of the enzyme molecules to the substrate are now deemed more realistic and reasonable since the conformations of all four enzyme molecules attained a stable state in the simulation that is closer to their natural conformations.
(A)
NIR binding state before (receptor: red, ligand: warm pink) and after (receptor: blue, ligand: purple) MDS with NO
2-
.
(B)
HAO binding state before (receptor: red, ligand: warm pink) and after (receptor: blue, ligand: purple) MDS with NH
2OH.
Figure 26.
Comparison of the enzyme's conformation with the relevant substrate docking data before and after simulation.
we reexamined the substrate's interaction with the amino acid residues in the enzyme molecule's pocket in the new binding state.
According to the analysis results, a nitrogen atom in the nitrite ion formed four hydrogen bonds with three amino acid residues in the simulated NAR pocket, two of which were contributed by the lysine at position 187, in addition to five amino acid residues that interacted hydrophobically with the nitrite ion from different angles
.
When compared to the docking results of the pre-simulation conformation, the number of hydrogen bonds formed and the total number of interacting amino acid residues between the enzyme conformation and the substrate were increased after the simulation, resulting in greater spatial site resistance to the substrate molecule
(Figure 28 A, B)
.
(A)
NIR amino acid residues react with NO
2-
NIR amino acid residues react with in three dimensions.
(B)
NIR amino acid residues react with NO
2-
in two dimensions.
Figure 27.
Docking results of NIR and NO
2-
after MDS.
The number of hydrogen bonds formed between amino acid residues and hydroxylamine in the simulated HAO pocket, as well as the number of amino acid residues with which there are interactions, are reduced when compared to the pre-simulation conformation, which may be related to the increased binding energy caused by molecular docking.
In the HAO pocket, the hydroxylamine molecule established four hydrogen bonds with three amino acid residues and hydrophobic interactions with threonine at position 63 and methionine at position 64.
.
(A)
HAO amino acid residues react with NH
2OH
in three dimensions.
(B)
HAO amino acid residues react with NH
2OH
in two dimensions.
Figure 28.
Docking results of HAO and NH
2OH after MDS.
Discussion
Through the above analysis, we
retrieved high-quality structure data
of the four enzyme molecules chosen for the subject for which no structural information is presently available, and
projected the placements of these enzyme molecules' pockets
.
Based on this positional information, we
obtained the details of the interactions
between the enzyme molecule and
the substrate through molecular docking and obtained the amino acid residues in the enzyme molecule that interact with the substrate
,
which is believed to play a key role in the catalytic process of its corresponding substrate ions. Therefore, the subsequent structural optimization and modification of this enzyme molecule should avoid the destruction of these active sites, which may lead to the enzyme loss of enzyme activity.
Furthermore, the study analyzed the flexible regions of the four enzyme molecules using molecular dynamics simulation, and these regions contributed significantly to the instability of the enzyme molecules, particularly the loop structure located in the middle part of the enzyme molecule, which may affect the stability of the conformation of the active center of the enzyme, and thus the subsequent modification
.
