In this model we describe our small biological switch inside Lactobacillus planetarium that only flips ON when two specific signals team up like a safety system that needs both a key and a code to open.
So , we are exploring a gene circuit designed as an AND gate to control our output which is CO-BERA mRNA. This gate is activated only when the environment becomes acidic (Ph below 6.9) and hydrogen peroxide (H2O2) is about 8q (MM) or more . Those double conditions are specific for our population which are people with severe uncontrolled asthma.
The activation of the first step starts when the pH equals 6.9; the pH begins to dip more, mimicking the environment in severe uncontrolled asthmatic patients . This drop activates P170_CP25 promoter to initiate LacR mRNA transcription that is translated to form LacR repressor as our first Hero. The guide Here is a pH sensitive promoter that boosts LacR mRNA transcription when pH value goes below 6.9.
Once LacR mRNA protein builds up, it shuts down the transcription of Rep repressor mRNA as LacR binds to its operator on P32 promoter that has the rule to produce Rep repressor. Normally, Rep repressor keeps CO-BERA locked away by repressing its transcription from pKat promoter. But when Rep repressor is out of the picture ,by the aid of Lac repressor, that lock opens. But that is only halfway.
The second key is the H2O2. As in normal pKat promoter is locked by a second lock which is per repressor, that has metal particles most of the time is Fe+2. This per repressor has a threshold in which it can sense H2O2, above this threshold H2O2 can oxidize the Fe+2 to be Fe+3 in this time the second lock is open and pKat can start transcription of our goal that is CO-BERA.
shows the double conditions that are essential to express CO-BERA are present so there is expression.
shows the double conditions that are essential to express CO-BERA are not present so there is no expression.
This system shows how synthetic biology can engineer bacteria able to sense and respond to dual signals, potentially for target therapies.
This Model boils down to 7 key equations describes the changes in our 7 main variables which are:
1- pH value
2- pH activation function for LacR induction
3- LacR mRNA
4- LacR protein
5- Rep repressor mRNA
6- Rep repressor protein
7- CO-BERA
This function models progressive acidification of the lung environment in asthma. pH evolves over time, starting at 6.9 and dropping gradually to simulate the environment in the lung in the condition of asthma.
This models our pH sensor which shows sigmoidal activation of LacR transcription as pH drops below 6.9. It is like a dimmer switch at normal pH lung environment it equals zero, but it gets the transcription done when the pH environment turns more acidic mimicking the asthma condition.
simulates the decrease in acidity mimicking the environment in the lung in severe uncontrolled asthma and its effect on F(pH) value over time.
This equation models the change for LacR mRNA (M), our starting point in this cascade:
-It increases due to basal transcription rate (V₀) plus the pH triggered boost (Vᵢ · F) then scaled by concentration of LacR gene in the cells(GR).
-It decreases due to natural degradation rate (dM) , cell dilution(μ) and translation rate to LacR repressor (C1), mimicking normal cellular clear up.
In this equations we see LacR repressor (LR) building up:
-It increases due to the translation rate for LacR protein (C1).
-It decreases due to natural degradation rate ( β).
This equation handles the Rep repressor mRNA (N):
-It increases due to the maximum transcription rate (V₂) and in presence of LacR repressor it inhibits the transcription by rate (Rₘₐₓ) scaled by sensitivity (S), then multiply by concentration of Rep repressor gene in the cells (G). In absence of LacR repressor we cap it at zero to avoid negative nonsense.
-It decreases due to natural degradation rate (dN) , cell dilution(μ) and translation rate to Rep repressor (C2) ,mimicking normal cellular clear up.
This equation shows Rep Repressor (P) building up:
-It increases due to the translation rate for Rep repressor (C₂).
-It decreases due to natural degradation rate ( γ).
Finally, our herooooooooooo CO-BERA :
CO-BERA production kicks in at rate (Tₘₐₓ) but that happens in high amounts of H₂O₂ above the promoter threshold that the halfway, the second half is the absence of Rep repressor then scaled by concentration of LacR gene in the cells(G) for consistency.
LacR active zone (Figure ): LacR transcription active zone is when pH is below 6.9 .
LacR expression (Figure ) simulates the response to the decrease in acidity as pH decreases below 6.9 it activates transcription of LacR mRNA that in turn is translated to LacR repressor.
LacR expression (Figure ) simulates the response to the increase in LacR repressor , that blocks the transcription of Rep mRNA that diminishes the Rep repressor as there is no production but only degradation.
AND Gate ON state (): AND Gate on when ( pH < 6.9) & (H2O2 > 8 micromole)
1.pH < 6.9, H₂O₂ < 8 μM: Minimal CO-BERA expression (~0.1 mM)
2.pH > 6.9, H₂O₂ Low CO-BERA expression (~2.3 mM)
3.pH < 6.9, H₂O₂ > 8 μM: Maximum CO-BERA expression (~22.3 mM)
4.pH > 6.9, H₂O₂ < 8 μM:Background expression only (~0.05 mM)
ON/OFF Ratio:446:1 (22.3 mM / 0.05 mM)
Basal Leakage Reduction: 95% compared to single-input controls
Response Time: 2-4 hours for full activation
Specificity:>10-fold selectivity for severe asthma conditions
At first project conditioning of CO-BERA expression was built on H2O2 dependant condition only but throw CO-LAB simulation of the H2O2 equations found that there is basal leakage that has a risk on the population as shown in the next figure.
Basal activity of pKatA (): Shows high basal expression of pKAT in absence of H2O2.
Basal activity of pKatA 2 (): Shows Basal activity at different H2O2 values.
After noticing this high basal activity, this credit goes to Model. We change to use our AND Gate conditioning system and the basal activity is nearly zero.
1.Dual Safety Lock:AND gate prevents activation from single stimuli
2.Threshold Calibration:pH and H₂O₂ thresholds set above normal lung variations
3.Response Kinetics:2-4 hour delay provides therapeutic window for intervention
4.Output Scaling:22.3 mM CO-BERA sufficient for downstream siRNA processing
Dual-fluorescent reporter system for independent pathway monitoring
Time-course analysis under simulated asthma conditions
Single vs. dual input specificity testing
pH/H₂O₂ threshold fine-tuning with clinical lung samples
1.Temporal Filtering: Sustained signal required (>2 hours) prevents transient activation
2.Reversible Logic:System returns to OFF when either input normalizes.
The final CO-BERA level when the AND Gate output is on = 22.321 (mM)
Within the validation of H2O2 key alone, we discovered that pKat promoter has basal leakage. To overcome this leakage we used another key, which is ph sensitive promoter, as a complementary condition to strengthen the lock and that decreased the basal activity to be near to zero.
Moreover, in some stress conditions H2O2 may increase to the threshold of its lock, which is not specific to our targets. On another hand after applying the complementary condition it helped us more specifically for asthmatic patients.
The validation of our loading system ,through molecular dynamics, was divided into two parts, combining an RNA-binding protein with an engineered RNA motif to ensure that our therapeutic RNA (CO-BERA) is efficiently packaged into bacterial membrane vesicles (BMVs). Our simulation consists of three models: a model for the Transmembrane proteins with the RNA-binding protein, a model for the CO-BERA, and finally a model for the endosomal escape protein.
The first part is the binding of the RNA-binding protein to the Lactobacillus membrane. This relies on the ribosomal protein L7Ae (RNA-binding protein) that identifies the K-turn motif C/D Box. We fused L7ae to a transmembrane protein to bind L7ae to bacterial membrane vesicles. We chose two different transmembrane proteins to compare, Foldase PrsA and DUF4811, using CHARMM-GUI to build the Lactobacillus plantarum membrane and GROMACS for molecular dynamics simulations.
The second part in our treatment is mainly dependent on small interfering RNA (siRNA) which will inhibit TSLP mRNA that will prevent it from being translated. This siRNA is in the CO-BERA scaffold which has a C/D Box sequence in its three terminals that will make it able to be bound to L7Ae (RNA-binding protein). L7Ae binds it to the linker protein and then to the transmembrane protein. This will deliver the whole complex to the Lactobacillus membrane to make the membrane vesicle and then be delivered to the lung cell. Therefore,we have docked all this complex in one pdb file, then we uploaded the pdb file on CHARMM-GUI to prepare the file for 5 nano second simulation on using OPENMM and AMBER tools simulated on google colab pro to detect the stability of the complex.
After viewing our loading system, we turned to validate our endosomal escape protein, Listeriolysin O, which is pH sensitive protein. Thus, we modeled two different pH levels, one acidic and the other neutral, and integrated them to the phagosomal membrane to test pore formation activity and stability.
After viewing our loading system, we turned to validate our endosomal escape protein, Listeriolysin O, which is pH sensitive protein. Thus, we modeled two forms of Listeriolysin O, one in Ph 6.5 and the other with a mutation in AminoAcid 461 to be Threonine instead of Leucine to be activated in the PH of the early endosome which is 6.5. This model will validate the activation of this PH level in the mutated form and the inactivation of the original form of the endosomal escape protein.
After determining the A.A sequences of all parts, it will be used in the SWISS model prediction. This computational process identifies the tertiary structure of a template protein using the Global Model Quality Estimate(GMQE) score during template search and selection. The higher the GMQE score, the more likely the correct structure, along with the coverage and the sequence identity percentage. The outcome of this step is the preparation of the two complexes, which will then be further processed by the CHARMM-GUI package to incorporate all simulation parameters and extract the 3D structure model of each complex.
The two complexes were uploaded to CHARMM-GUI Membrane Builder. Then we selected the full protein chain and oriented it automatically using the PPM 2.0 tool, which places transmembrane proteins correctly within the bilayer.
To mimic the natural membrane of L. plantarum, we used experimental lipid composition data from Pierart et al. (2015. The bilayer was built with the following study in the references.
The membrane–protein system was solvated with water layers (22.5 Å thickness), and physiological ion concentration (0.15 M KCl) was added for stability and charge neutralization.
We used the CHARMM36m force field, which provides accurate parameters for proteins, lipids, and ions. Hydrogen mass repartitioning was enabled to allow longer time steps in molecular dynamics simulations.
● Step 5: Final System
CHARMM-GUI generated the full input files for the GROMACS simulation. The final box size was ~111 × 111 × 134 Å, containing:
- Either (DUF4811–Linker–L7Ae) or (FOLDASE-Linker-L7Ae)
- ~20,000 water molecules
- K⁺ and Cl⁻ ions
After finishing these steps in both complexes.
We are ready now to extract and interpret the RMSD & RMSF graphs of each complex to compare them. The graphs were extracted by Gromacs, and then VMD will be used for visualizing their motion.
1- DUF4811 -Linker - RNA binding protein L7ae complex simulation got the best simulation results as it showed higher stability and rigidity compared to Foldase PrsA -Linker - RNA binding protein L7ae complex. When simulated at 5ns, the root mean square deviation (RMSD) ranged from 1.0 Å to 2.0 Å, which is efficiently stable within the membrane. Also, the Root Mean Squared Fluctuation (RMSF) ranged from 0.1 Å to 0.2 Å, which ensured complex inflexibility.
2-Foldase PrsA complex simulation showed a higher result than DUF4811 complex, as when simulated at 5ns, the RMSD ranged from 5.0 Å to 9.0 Å, which is a very high range and shows that the complex is unstable. RMSF ranged from 0.2 Å to 0.5 Å, and it showed a high increase in the beginning, which validates that the complex is fluctuant and flexible.
Finally, this RMSD comparison shows the huge variability between the two transmembrane proteins, as the DUF4811 complex expressed a very narrow peak, low structural variability, and high stability compared to the Foldase PrsA complex, which has a broader distribution starting from 2 Å to 10 Å, much higher structural variability, and flexibility.
Comprehensive Analysis of TMP2 Transmembrane Integration
Detailed Analysis of TMP3 Membrane Integration Challenges
Integrated Assessment of All Three L7AE Fusion Complexes
To prepare each complex, we began with the Amino Acid sequence of each protein.
After determining the A.A sequences of all parts, it will be used in the SWISS model prediction. This computational process identifies the tertiary structure of a template protein using the Global Model Quality Estimate(GMQE) score during template search and selection. The higher the GMQE score, the more likely the correct structure, along with the coverage and the sequence identity percentage. The outcome of this step is the preparation of the two complexes, which will then be further processed by the CHARMM-GUI package to incorporate all simulation parameters and extract the 3D structure model of each complex.
After obtaining the two complex structures and coordinates (.pdb files) from SWISS MODEL, we turned to CHARMM-GUI to prepare them with our parameters. It underwent several steps to achieve these points, including placing our engineered protein inside a Built L. plantarum membrane, testing lipid interactions, protein stability, and orientation, and preparing the system for long molecular dynamics simulations. These simulations will help us evaluate protein stability and linker flexibility. The upcoming steps illustrate how we built the membrane and prepared the files for the simulation.
● Step 1: Importing the 3D structures (.pdb files) into CHARMM-GUI
The two complexes were uploaded to CHARMM-GUI Membrane Builder. Then we selected the full protein chain and oriented it automatically using the PPM 2.0 tool, which places transmembrane proteins correctly within the bilayer.
● Step 2: Membrane Composition
To mimic the natural membrane of L. plantarum, we used experimental lipid composition data from Pierart et al. (2015. The bilayer was built with the following study in the references.
● Step 3: Solvation and Ions importing
The membrane–protein system was solvated with water layers (22.5 Å thickness), and physiological ion concentration (0.15 M KCl) was added for stability and charge neutralization.
● Step 4: Force Field
We used the CHARMM36m force field, which provides accurate parameters for proteins, lipids, and ions. Hydrogen mass repartitioning was enabled to allow longer time steps in molecular dynamics simulations.
● Step 5: Final System
CHARMM-GUI generated the full input files for the GROMACS simulation. The final box size was ~111 × 111 × 134 Å, containing:
- Either (DUF4811–Linker–L7Ae) or (FOLDASE-Linker-L7Ae)
- ~20,000 water molecules
- K⁺ and Cl⁻ ions
When files are prepared by CHARMM-GUI, we are now ready to simulate our two complexes using the GROMACS 2025.2 version.
The molecular dynamics simulations were performed using GROMACS 2025.1 on an Ubuntu Linux system. This phase is critical for generating the output trajectory files, which are essential for conducting stability analysis, specifically showing the Root Mean Square Deviation (RMSD) and Root Mean Square Fluctuation (RMSF) graphs. These metrics are key to comparing the stability of the two L7AE-linker-transmembrane protein interactions. The simulation workflow began by utilizing the prepared GROMACS input files from CHARMM-GUI. These inputs include: parameter files (.mdp) for energy minimization, equilibration, and production; topology files (.top); an index file (.ndx); and the initial coordinate file (.gro).
After installing GROMACS on Ubuntu, we performed the simulation steps as following: one minimization step to prevent simulation instability and to remove any unnatural bond lengths (steric clashes); six distinct equilibration steps were conducted to gradually bring the system to the target temperature (NVT), achieve correct density and pressure (NPT); and one production step to observe, analyze the dynamic process. The final trajectory file from the production run is the key output used for conducting the stability analysis. By repeating this entire simulation procedure for the L7AE RNA-binding protein linked to each of the two different target proteins, our succeeding analysis revealed that the complex containing (L7AE-linker-DUF4811) exhibited the greatest stability.
● The first step of the simulation: the minimization step
gmx grompp -f step6.0_minimization.mdp -c step5_input.gro -p topol.top -o em.tpr -r step5_input.gro -n index.ndx
gmx mdrun -v -deffnm em
● Second step: the 6 equilibration steps to prepare the complex environment
gmx grompp -f step6.1_minimization.mdp -c em.gro -p topol.top -o nvt.tpr -r em.gro -n index.ndx
gmx mdrun -v -deffnm nvt
gmx grompp -f step6.2_equilibration.mdp -c nvt.gro -p topol.top -o npt2.tpr -r nvt.gro -n index.ndx
gmx mdrun -v -deffnm npt2
gmx grompp -f step6.3_equilibration.mdp -c npt2.gro -p topol.top -o npt3.tpr -r npt2.gro -n index.ndx
gmx mdrun -v -deffnm npt3
gmx grompp -f step6.4_equilibration.mdp -c npt3.gro -p topol.top -o npt4.tpr -r npt3.gro -n index.ndx
gmx mdrun -v -deffnm npt4
gmx grompp -f step6.5_equilibration.mdp -c npt4.gro -p topol.top -o npt5.tpr -r npt4.gro -n index.ndx
gmx mdrun -v -deffnm npt5
gmx grompp -f step6.6_equilibration.mdp -c npt5.gro -p topol.top -o npt6.tpr -r npt5.gro -n index.ndx
gmx mdrun -v -deffnm npt6
● Finally, the production step:
gmx grompp -f step7_production.mdp -c npt6.gro -p topol.top -o md.tpr -n index.ndx
gmx mdrun -v -deffnm md
After finishing these steps in both complexes.
We are ready now to extract and interpret the RMSD & RMSF graphs of each complex to compare them. The graphs were extracted by Gromacs, and then VMD will be used for visualizing their motion.
1- DUF4811 -Linker - RNA binding protein L7ae complex simulation got the best simulation results as it showed higher stability and rigidity compared to Foldase PrsA -Linker - RNA binding protein L7ae complex. When simulated at 5ns, the root mean square deviation (RMSD) ranged from 1.0 Å to 2.0 Å, which is efficiently stable within the membrane. Also, the Root Mean Squared Fluctuation (RMSF) ranged from 0.1 Å to 0.2 Å, which ensured complex inflexibility.
2-Foldase PrsA complex simulation showed a higher result than DUF4811 complex, as when simulated at 5ns, the RMSD ranged from 5.0 Å to 9.0 Å, which is a very high range and shows that the complex is unstable. RMSF ranged from 0.2 Å to 0.5 Å, and it showed a high increase in the beginning, which validates that the complex is fluctuant and flexible.
Finally, this RMSD comparison shows the huge variability between the two transmembrane proteins, as the DUF4811 complex expressed a very narrow peak, low structural variability, and high stability compared to the Foldase PrsA complex, which has a broader distribution starting from 2 Å to 10 Å, much higher structural variability, and flexibility.
This model shows how we tested the whole complex of siRNA in the CO-BERA scaffold with the protein complex that will deliver it to the membrane of the bacteria to make membrane vesicle and then be delivered to the lung cells.
To bind our CO-BERA sequence with the L7Ae (RNA-binding protein) we need C/D box sequence, to determine how many copies of the C/D Box we need; so we compared between 4 different C/D box copies from 1 copy to 4 copies.
First we upload the fasta sequence of each one on RNA FOLD server and compare between them according to Ensemble Free Energy (average free energy of all possible RNA secondary structures), Frequency of MFE Structure (percentage of RNA molecules would actually adopt the predicted minimum free energy),and Ensemble Diversity (how many different structures the RNA can adopt) and gave us this results :
1.One copy of C/D Box:
fig() this figure shows the CO-BERA with one C/D Box copy visualization and its mountain plot
The free energy of the thermodynamic ensemble is -130.86 kcal/mol.
The frequency of the MFE structure in the ensemble is 0.01 %.
The ensemble diversity is 45.57
2.Two copies of C/D Box:
fig() this figure shows the CO-BERA with two C/D Box copies visualization and its mountain plot
The free energy of the thermodynamic ensemble is -141.87 kcal/mol.
The frequency of the MFE structure in the ensemble is 0.00 %.
The ensemble diversity is 47.72 .
3.Three copiesof C/D Box:
fig() this figure shows the CO-BERA with three C/D Box copies visualization and its mountain plot
The free energy of the thermodynamic ensemble is -158.93 kcal/mol.
The frequency of the MFE structure in the ensemble is 0.00 %.
The ensemble diversity is 54.71
4.Four copiesof C/D Box:
fig() this figure shows the CO-BERA with four C/D Box copies visualization and its mountain plot
The free energy of the thermodynamic ensemble is -170.41 kcal/mol.
The frequency of the MFE structure in the ensemble is 0.00 %.
The ensemble diversity is 88.04 .
fig() Bar chart that shows the values of ensemble diversity of the four different copies.
After interpreting the results 1 and 2 copies of C/D Box have relatively lower ensemble diversity than others, giving us a more predictable structure that RNA can adopt.
To determine which one would be the choice we need to dock them with L7Ae to compare the docking score of each one, this had been done using HDOCK server and gave us the following results:
1.one copy of C/D Box:
fig() Gif shows visualization of CO-BERA with one C/D box copy docked with L7Ae.
Docking score: -284.05
2.two copies of C/D Box
fig() Gif shows visualization of CO-BERA with two C/D box copies docked with L7Ae.
Docking score: -234.54
1.Three copies of C/D Box:
fig() Gif shows visualization of CO-BERA with three C/D box copies docked with L7Ae.
Docking score: -268.95
2.three copies of C/D Box
fig() Gif shows visualization of CO-BERA with four C/D box copies docked with L7Ae.
Docking score:-249.59
fig() Bar chart shows the different docking scores of the four complexes
Using HDOCK revealed that the docking score of 1 copy of C/D box is significantly higher than 2 copies, which made us choose 1 copy of C/D box as the most stable one for our simulation.
We made our complex by using Alphafold server by inserting the sequence of the CO-BERA-CDbox complex with the complex protein which consists of RNA-binding protein(l7ae)- linker protein-transmembrane protein (DUF4811 protein) as this form the docked file in CIF form which have been changed to pdf form using science codon server.
This formed the docked file in CIF form which has been changed to pdf form using science codon server.
As in the previous model the pdb file has been uploaded to CHARMM GUI server to prepare the pdb file for the simulation by which:
CHARMM-GUI detects: 3 protein chains (for L7Ae ,linker and DUF4811) and RNA (for CO-BERA and C/D Box)
We have selected KCLions with 0.15 concentration to neutralize the system and achieve the desired system concentration
In this step we solved our system in water and PH mimic the physiological condition. We have chosen an octahedral box.
We applied for Proteins → ff14SB ,RNA → OL3 Water → TIP3P and as an input generation we have chosen AMBER and GROMACS
Minimization input : amber/step4.0_minimzation.mdin, equilibration input: amber/step4.1_equilibration.mdin and dynamics input amber/step5_production.mdin
After that we ran the output files on google colab pro template that run OPENMM engine + AMEBR tools this will help us to detect the stability of our complex by doing Root Mean Square Deviation (RMSD) graph and Root Mean Square Fluctuation (RMSF) graph these analysis verify that the simulation captures realistic binding interaction, the simulation work flow began by utilizing the prepared AMBER input form CHARMM-GUI.
After installing the OPENMM engine and AMBER tools on GOOGLE COLAB PRO, we performed the simulation steps as the previous model beginning from minimization. We make the system undergo three steps of minimization with gradually stricter tolerances (gentle → moderate → final). Equilibration, the system is heated gradually from 0 → 100K → 300K to bring the system to the target temperature (NVT), achieve correct density and pressure (NPT); and the production run which produces the final trajectory file which is used for conducting the stability analysis.
● The first step: minimization step:
illu
● The second step: equilibration step:
To mimic the natural membrane of L. plantarum, we used experimental lipid composition data from Pierart et al. (2015. The bilayer was built with the following study in the references.
● The third step production step:
illu
After finishing these steps in both complexes.
We are ready now to extract and interpret the RMSD & RMSF graphs of each complex to compare them
The graph shows stabilized Root Mean Square Deviation (RMSD) stabilized between 3.0-5.0 Å (0.30-0.50 nm), which is within the acceptable range for siRNA-protein complex, and it falls in the optimal range for the functionality of our complex as the duplex siRNA will unwind forming 2 siRNA strand, also is essential for RISC loading and Dicer processing.
The RMSF graph show values variation from 0.8-5.5 Å, the RMSF values of 0.8-2.0 Å representing stable duplex core segments that maintain overall structural integrity during processing, RMSF values of 2.0-3.0 Å indicating controlled duplex breathing necessary for thermodynamic asymmetry assessment the extremely flexibility peak at residue 50 (RMSF=5.5 Å) which corresponds to the duplex cleavage site by DICER, the higher flexibility the higher enzyme accessibility and strand displacement. So, the RMSF balanced flexibility shows effective RNA interference pathway activation
First we identified the listeriolysin O Protein ID(4CDB) from the RSCB protein data bank. We went to build the membrane on the CHARMM-GUI Package
● Step 1: Importing the 3D structures (.pdb files) into CHARMM-GUI
Listeriolysin O ID was uploaded to CHARMM-GUI Membrane Builder. Then we selected the full protein chain and oriented it automatically using the PPM 2.0 tool, which places protein correctly within the membrane. We made this step twice, once with the protein in PH 6.5 ,and the other with mutation in AMINO ACID 461 to be threonine instead of leucine also in PH 6.5.
● Step 2: Membrane Composition
To mimic the natural membrane of a phagosome, we used experimental lipid composition data from Karin A. Zemski Berry, Robert C. Murphy, Beata Kosmider, and Robert J. Mason (2017. The bilayer was built with the Following study in the references.
● Step 3: Solvation and Ions importing
The membrane–protein system was solvated with water layers (22.5 Å thickness), and physiological ion concentration (0.15 M KCl) was added for stability and charge neutralization.
● Step 4: Force Field
We used the CHARMM36m force field, which provides accurate parameters for proteins, lipids, and ions. Hydrogen mass repartitioning was enabled to allow longer time steps in molecular dynamics simulations.
● Step 5: Final System
CHARMM-GUI generated the full input files for the GROMACS simulation. The final box size was ~111 × 111 × 134 Å, containing:
- Either Listeriolysin O in PH 6.5 or mutated Listeriolysin O in PH 6.5
- ~20,000 water molecules
- K⁺ and Cl⁻ ions
In this step, we used the Gromacs software to simulate the topology files downloaded from the CHARMM-GUI package. We used the same steps mentioned in - under the Gromacs tutorial. After finishing the simulation, we will analyze the trajectory files by using GROMACS and then the simulation visualization by PYMOL.
To prepare Listeriolysin O, we began with the AminoAcid sequence.
After determining the A.A sequences of all parts, it will be used in the SWISS model prediction. This computational process identifies the tertiary structure of a template protein using the Global Model Quality Estimate(GMQE) during template search and selection. The higher the GMQE score, the more likely the correct oligomer structure. The outcome of this step is the preparation of the protein structure and coordinates file to be simulated by the CHARMM-GUI package.
● Step 1: Importing the 3D structures (.pdb files) into CHARMM-GUI
The Listeriolysin O pdb file was uploaded to CHARMM-GUI Membrane Builder. Then we selected the full protein chain and oriented it automatically using the PPM 2.0 tool, which places the protein correctly within the bilayer. We made this step twice, once in PH 6 and the other in PH 7.4.
● Step 2: Membrane Composition
To mimic the natural membrane of a phagosome, we used experimental lipid composition data from Karin A. Zemski Berry, Robert C. Murphy, Beata Kosmider, and Robert J. Mason (2017. The bilayer was built with the Following study in the references.
● Step 3: Solvation and Ions importing
The membrane–protein system was solvated with water layers (22.5 Å thickness), and physiological ion concentration (0.15 M KCl) was added for stability and charge neutralization.
● Step 4: Force Field
We used the CHARMM36m force field, which provides accurate parameters for proteins, lipids, and ions. Hydrogen mass repartitioning was enabled to allow longer time steps in molecular dynamics simulations.
● Step 5: Final System
CHARMM-GUI generated the full input files for the GROMACS simulation. The final box size was ~111 × 111 × 134 Å, containing:
- Either Listeriolysin O in PH 6 or Listeriolysin O in PH 7.
- ~20,000 water molecules
- K⁺ and Cl⁻ ions
In this step, we used the Gromacs software to simulate the topology files downloaded from the CHARMM-GUI package. We used the same steps mentioned in MODEL 1 under the Gromacs tutorial. After finishing the simulation, we will analyze the trajectory files by using GROMACS and then the simulation visualization by PYMOL.
