Introduction

We simulated gene expression driven by P43, a constitutive promoter, using Ordinary differential Equations.

A model for constitutive gene expression was developed in MATLAB utilizing ordinary differential equations (ODEs). The core components of the system, namely mRNA and protein species, were defined to reflect the gene expression dynamics. The transcription process was implemented as a first-order reaction, wherein mRNA is produced at a constant transcription rate, characteristic of constitutive expression. Additionally, the degradation of both mRNA and protein was incorporated as first-order processes, each defined by their respective degradation rates. The system of ODEs was formulated to describe the dynamics of mRNA and protein concentrations over time. MATLAB's ode45 solver was employed to simulate the model across a specified time span, enabling the analysis of the resulting mRNA and protein concentrations. The outputs were plotted to visualize the behavior of the gene expression system.

The data for rate constants were derived from the results in the paper “Engineering an inducible gene expression system for Bacillus subtilis from a strong constitutive promoter and a theophylline-activated synthetic riboswitch” [4].

Fig: Diagram of the model for constitutive gene expression
Fig: The graph is obtained for proteins from constitutive gene expression from the model’s simulation. X-axis is time in hours and Y-axis is number of molecules.
Fig: Equations used for making the model.

Induced Gene Expression Simulation for PbacA

For PbacA, a repressible promoter, induced gene expression using ODE was simulated.

A model for induced gene expression was developed in MATLAB using ordinary differential equations (ODEs) to capture the dynamics of gene activation in response to an inducer. The system comprised key components, including mRNA and protein species, as well as the transcription factor (inducer) concentration. The transcription process was defined using a Hill equation to account for the cooperative binding of the transcription factor to the promoter, reflecting the dynamics of induced expression. Additionally, the degradation of mRNA and protein was included as first-order processes, each characterized by their respective degradation rates. A system of ODEs was formulated to describe the time-dependent behavior of mRNA and protein concentrations, influenced by the varying concentration of the inducer. The model was simulated over a specified time span using MATLAB's ode45 solver, allowing for the analysis of the resulting mRNA and protein levels in response to changes in the inducer concentration. The simulation outputs were plotted to visualize the dynamics of the gene expression system.

The data used for the model was derived from the results obtained in the paper - “Construction and Characterization of a Gradient Strength Promoter Library for Fine-Tuned Gene Expression in Bacillus licheniformis” [5].

Fig : Diagram of the induced gene expression model.
Fig: The graph obtained for proteins translated in induced gene expression from the model’s simulation. X-axis is time in hours and Y-axis is number of molecules.
Fig: Equations for the induced gene expression model

We intend to replace this natural promoter with the synthetic cumate induced promoter constructed by combining the strong constitutive Bacillus promoter Pveg with regulatory elements of the Pseudomonas putida, CymR repressor, and its operator sequence CuO [2].

For our implementation, cumate will be added into our formulation of wettable powder which will be taken up by our chassis passively. The formulation (chassis + Cumate) will be sprayed onto the coffee leaves in the plantation to distribute the chassis for its antifungal action against the HV. Once the Cumate has degraded or washed away the bacteria will no longer be able to survive in that environment. Cumate is a biodegradable inexpensive lab chemical. It is degraded by several microorganisms especially alphaproteobacterial like Rhodococcus,Pseudomonas, and Sphingomonas and is not toxic to the environment [3][4][5].

Docking

Promoter Strength and Sigma Factor Analysis

A promoter's strength is primarily classified by the stability of its open complex and the binding affinity to the RNAP. The binding affinity is measured through docking with the sigma factors.

Since the structure for Sigma B + RNAP (holoenzyme) was not available on any available databases, we superimposed the structure of Holoenzyme Sigma A onto Sigma B, whose structure was predicted using AlphaFold, without the conformational changes due to RNAP [3].

The structure for holoenzyme Sigma A was obtained from the PDB database (7CKQ). This was visualized using PyMOL alongside the RNA polymerase to assess the similarity between Sigma B and Sigma A. This analysis indicated whether the conformational change in the sigma factor was accounted for in the Sigma B obtained from the PDB database.

Core Promoter Visualization and Docking Analysis

Core promoter regions of P43 optimized, P43 unoptimized, PbacA optimized, and PbacA unoptimized were visualized using PyMOL. In Bacillus subtilis, DNA is in the A form during the sporulation stage; hence, PbacA fragments were modeled as A helices, while P43 fragments were modeled as B helices [1][2].

Based on the superimposition of Sigma A (272 residues) onto Sigma B (262 residues), 233 residues were structurally aligned. Due to the high similarity between the two structures, we assume the two structures to be in a similar conformational state.

Docking of sigma factors to the core regions of the promoters P43 (both unoptimized and optimized) and PbacA (both unoptimized and optimized) at the –10 and –35 regions was performed to evaluate whether the optimization exhibited improved binding affinity towards the sigma factors, which would imply stronger gene expression.

The optimization performed involved changing the –35 and –10 regions to the predicted consensus sequence.

P43 unoptimized and optimized core region fragments were docked with Sigma Factor A from its RNAP holoenzyme (PDB ID: 7CKQ). The same two fragments were also docked with Sigma Factor B, whose structure was obtained from AlphaFold (AF-A0A6M4JFM2-F1-v4). PbacA optimized and unoptimized core regions were docked with Sigma Factor B.

Fig. CymR stays bound to the operator repressing the promoter (-41:+8). Cumate lifts the bound repressor in the presence of which galE is expressed
Fig: difference between PbacA optimised and pbacA unoptimised

Docking Scores and Binding Affinities

The docking score directly correlates with binding affinity. More negative Gibbs' free energy implies stronger binding. From the docking scores, we can see that:

  • P43 unoptimized showed better binding affinity compared to P43 optimized with Sigma A and Sigma B.
  • P43 optimized shows more binding affinity towards Sigma B.
  • PbacA optimized core regions show better binding affinity towards Sigma Factor B.
  • The binding affinity of PbacA optimized was better than that of P43 unoptimized and P43 optimized, while PbacA unoptimized had a lesser binding affinity compared to P43 unoptimized with Sigma Factor B.
P43 Optimised + Sigma A
Rank 1 2 3 4 5 6 7 8 9 10
Docking Score -204.52 -196.27 -195.33 -189.22 -188.04 -187.96 -186.45 -186.07 -183.36 -182.8
Confidence Score 0.7485 0.7161 0.7123 0.6866 0.6815 0.6812 0.6746 0.6729 0.6609 0.6584
Ligand rmsd (Å) 301.14 253 347.4 274.38 289.19 335.12 284.53 324.15 254.53 308.5

P43 Unoptimised + Sigma A
Rank 1 2 3 4 5 6 7 8 9 10
Docking Score -243.3 -197.63 -194.77 -191.9 -188.94 -188.9 -186.9 -185.01 -184.58 -183.74
Confidence Score 0.866 0.7216 0.71 0.698 0.6854 0.6852 0.6766 0.6682 0.6663 0.6626
Ligand rmsd (Å) 254.45 304.65 333.12 270.16 264.77 247.3 277.32 272.69 300.31 315.96

P43 Unoptimised + Sigma B
Rank 1 2 3 4 5 6 7 8 9 10
Docking Score -240.37 -239.03 -213.4 -204.19 -202.98 -201.88 -201.33 -199.87 -195 -192.06
Confidence Score 0.859 0.8558 0.7804 0.7472 0.7426 0.7384 0.7363 0.7305 0.7109 0.6987
Ligand rmsd (Å) 84.62 81.51 55.12 83.81 114.25 76.66 114.28 71.45 109.95 91.42

P43 Optimised + Sigma B
Rank 1 2 3 4 5 6 7 8 9 10
Docking Score -235.5 -234.24 -214.96 -213.53 -212.63 -204.04 -203.81 -203.77 -199.03 -199.01
Confidence Score 0.8468 0.8435 0.7857 0.7808 0.7777 0.7466 0.7458 0.7456 0.7272 0.7271
Ligand rmsd (Å) 113.89 101.09 120.59 52.72 107.89 116.65 110.12 93.94 137.24 82.15

PbacA Unoptimised + Sigma B
Rank 1 2 3 4 5 6 7 8 9 10
Docking Score -229.9 -227.85 -227.35 -225.14 -223.5 -217.98 -212.95 -210.68 -210.45 -208.53
Confidence Score 0.8317 0.8259 0.8245 0.818 0.8131 0.7957 0.7789 0.7709 0.7701 0.7633
Ligand rmsd (Å) 65.8 80.98 68.21 71.54 72.25 79.69 71.46 55.29 72.97 81.29

PbacA Optimised + Sigma B
Rank 1 2 3 4 5 6 7 8 9 10
Docking Score -245.43 -225.89 -207.36 -204.87 -203.11 -200.97 -200.26 -199.02 -198.85 -197.83
Confidence Score 0.8709 0.8202 0.759 0.7498 0.7431 0.7349 0.7321 0.7272 0.7265 0.7224
Ligand rmsd (Å) 232.69 236.55 267.44 244.18 219.93 252.69 266.93 254.95 212.36 237.3

Conclusion: From the results, we concluded that the optimization for P43 is not effective, while the optimization of PbacA was effective and perhaps necessary for the functioning of our dual promoter.

Molecular Dynamic Studies

Molecular Dynamics (MD) Simulation was carried out for the Protein-DNA complexes using the Gromacs v2024.1 software with AMBER99 forcefield for 100 ns time.

Results

P43 optimised-SigA

Fig: Molecular Dock of P43 Optimised with Sigma A
Fig: RMSD Values for P43 Optimised with Sigma A

P43 optimised-SigB

Fig: Molecular Dock of P43 Optimised with Sigma B
Fig: RMSD Values for P43 Optimised with Sigma B

P43 unoptimised-SigA

Fig: Molecular Dock of P43 Unoptimised with Sigma A
Fig: RMSD Values for P43 Unoptimised with Sigma A

P43 unoptimised-SigB

Fig: Molecular Dock of P43 Unoptimised with Sigma B
Fig: RMSD Values for P43 Unoptimised with Sigma B

PbacA optimised-SigB

Fig: Molecular Dock of PbacA Optimised with Sigma B
Fig: RMSD Values for PbacA Optimised with Sigma B

PbacA unoptimised-SigB

Fig: Molecular Dock of PbacA Unoptimised with Sigma B
Fig: RMSD Values for PbacA Unoptimised with Sigma B

Interpretation

The Root Mean Square Deviation (RMSD) of the Protein-DNA complexes illustrate that the optimized DNA-Protein complexes are more stable than the unoptimized complexes. All the RMSD graphs are deviating within 2 Angstroms, which shows that all the complexes are stable. However, optimized complexes are more stable than the unoptimized complexes.

Conclusion

Optimized DNA exhibit better interactions with the protein compared to the unoptimized DNA.

References

  1. Setlow, P. (1992). DNA in dormant spores of Bacillus species is in an A-like conformation. Molecular Microbiology, 6(5), 563–567. https://doi.org/10.1111/j.1365-2958.1992.tb01501.x
  2. Jameson, K. H., & Wilkinson, A. J. (2017). Control of Initiation of DNA Replication in Bacillus subtilis and Escherichia coli. Genes, 8(1). https://doi.org/10.3390/genes8010022
  3. Fang, C.L., Zhang, Y. The cryo-EM structure of B. subtilis BmrR transcription activation complex PDB DOI: https://doi.org/10.2210/pdb7CKQ/pdb
  4. Cui, W., Han, L., Cheng, J., Liu, Z., Zhou, L., Guo, J., & Zhou, Z. (2016). Engineering an inducible gene expression system for Bacillus subtilis from a strong constitutive promoter and a theophylline-activated synthetic riboswitch. Microbial Cell Factories, 15(1). https://doi.org/10.1186/s12934-016-0599-z
  5. Rao, Y., Li, P., Xie, X., Li, J., Liao, Y., Ma, X., Cai, D., & Chen, S. (2021). Construction and Characterization of a Gradient Strength Promoter Library for Fine-Tuned Gene Expression in Bacillus licheniformis. ACS Synthetic Biology, 10(9), 2331–2339. https://doi.org/10.1021/acssynbio.1c00242