The goal of the modelling section was to test the predicted behaviour of the gene circuit and feedback into its design while testing its behaviour from different perspectives. Our modelling comprised three main sections. To clarify, the system is made up of the shinogen, the shinorine-producing cluster of genes extracted from Anabaena variabilis and the thanogen, the novel killswitch designed to take inputs of glucose and green light – a proxy for UV. These sections included the molecular modelling and the dynamics of the gene circuit within a single bacterium; biosafety modelling, considering the evolutionary stability of the thanogen to genetic drift and selection; and the macroscopic modelling, testing the behaviour of an entire transformed bacterial population carrying the gene circuit, as is the intended environment.
In the first, the molecular model of the behaviour of the gene circuit was dissected, taking account of the designed properties we sought in the design and research section of our engineering cycle and the intrinsic properties of the system (diffusion rates, association/dissociation, promoter dynamics). Therefore, the ultimate target this section was to paint a picture of the behaviour of the gene circuit overall and feed into the other two sections. Systems of ODEs were built for the thanogen and shinogen independently and deterministic and stochastic models produced using kinetic parameters derived from the literature and approximated where unavailable.
For the biosafety modelling, we designed a novel mathematical model to test the evolutionary stability of our thanogen. This involved initialising the system of ODEs derived from the microscopic model the thanogen and inputting its kinetic parameters as a function of the individual parts that made up its circuit – i.e. promoters, CDS and degradation tags. Using the length (in base-pairs) of each part in the gene circuit, we could approximate the likelihood of a mutation occurring in a given portion of the gene circuit and thus test which parts in the thanogen were most susceptible to dysfunction across evolutionary time. We then output this as a function of number of generations.
The molecular modelling discussed in the first section is very useful for informing gene circuit design and optimisation. However, this provides limited insight into how our technology might perform in any realistic application. Where very large numbers of bacteria are combined, simultaneous consideration of their individual behaviour soon becomes unfeasible. Our transformed E.coli bacteria will need to be contained within some gel or lotion formulation which can realistically be applied to the skin, in which the bacteria must survive and grow. The bacteria are likely to interact with each other, the environment, the ‘gel’, and with the other microflora contained in the microbiome of the skin. Intricate macroscopic behaviour soon arises from the simple microscopic rules and to make any predictions a statistical approach to the modelling process must be used. This section therefore looked into the behaviour of the system as a whole, taking account for bacterial proliferation and death while using the derived concentration of shinorine to model net protection to test shinescreen’s efficiency in UV protection.
Background – Biosafety Section
Mutation is an inevitable part of biological systems. Genetic mutations risk non-functionalisation of the gene network and a loss of stability in the population over time. The introduced gene circuit could confer a selective disadvantage to the carrier due to metabolic cost or gene toxicity and thus be selected against and favouring loss of the genetic material. This is further compounded by the instability of plasmid-borne gene where post-segregational loss of material is common during cell division. Therefore, evolutionary stability of the thanogen is crucial for both minimising the risk of undesirable long-term survival of the chassis bacterium outside the epidermis and uniform behaviour across the population.
Our product is a probiotic and is designed to be applied to the skin. Therefore, it is vitally important that our safety mechanisms are extremely robust. Consequences of our safety mechanisms failing could include:
- Permanent integration into the skin, or a permanent change in the composition of the skin microbiome
- Survival of the bacteria away from the skin and integration into the environment
During our interviews with various researchers, the issues of biosafety and avoidance of integration into the environment came up repeatedly. These interviews are summarised in our Human Practices page. Due to this input, we redesigned our plasmid to make it more robust and evolutionarily stable, and therefore safer. We also decided to study the evolutionarily stability of our kill switch through modelling.
To that end, we devised an evolutionary mutation algorithm we called the ‘Genetic Algorithm’ to test the resistance of the thanogen specifically to point mutations across n generations. Essentially, the algorithm uses the kinetic parameters from the ODEs of the thanogenic model, the net number of bps of each contributing part in the gene circuit and a belt of approximate mutation rates per generation of E. coli as inputs. A matrix of bacteria of size S is initialised with the starting parameters, as specified in the thanogen model. It is assumed that at t=0, ie when the bacteria are first applied to the skin, that the killswitch is operating optimally – the initial parameters implicit in the gene network are uniform across the population. Then, the Genetic Algorithm is cycled for n generations whereby seven essential processes occur:
- 1. A random number from a Poisson distribution with l = 1 is selected (see fig 1); this determines the number of bacteria within a population to undergo mutation. This is coupled with an isotropic random-number selector which, using the predicted number of bacteria to mutate, selects individual bacteria in the population by ID.
- 2. A second number from a Poisson distribution assigns the number of intracellular theoretical mutation sites to occur, adjusted to the number of bp of the gene circuit (total size of the plasmids times copy number) and a subsequent random-number generator selects specific loci in the circuit to mutate.
These loci are then fed into the ‘mutation machine’ which determines the kinetic parameter impacted by the mutation. Parts for which kinetic parameters were correlated were the promoters (transcription rates), CDS (binding and transcriptional activation/repression) and ssrA degradation tags (degradation). Depending on whether the part were a CDS or regulatory sequence, the the known probability that a given mutation site should be silent, missense, nonsense or a frameshift mutation could be translated into a mutational effect size.
- 4. The corresponding parameter is adjusted by summing the parameter at t-1 and the per nucleotide parameter contribution times a randomly generated normally distributed deviation.
- 5. After all mutation outcomes have been assigned to mutation sites for bacteria within the population at generation n, the bacteria undergo a cell division event and 50% of daughter cells either randomly (if drift simulation is being run) or non-randomly according to fitness are eliminated from the population. This changes t he distribution of alleles within the population.
- 6. The system of ODEs is re-integrated per bacterium or averaged across the population and the concentrations of A2B4 and endonuclease outputted by generation.
The former approach allows Ge, the generation at which the escape proportion =5% to be determined. The above cycle is then repeated.
- The killswitch is working at an optimum efficiency at generation 0
- The parameters in the ODE map directly on to different parts of the killswitch
- Daughter cells inherit all parameters from their parent
- The number of mutations which occur is proportional to the magnitude of the DNA
- Changes in part parameters are independent and epistatic interactions are negligible
- For the drift model, it is assumed that the killswitch is effectively silent during the cell cycle
The model was further designed to account for multiple replicates, accounting for the multiple copies of the thanogen.
The algorithm approach was used across our present thanogen. In minute population sizes <100 and a high mutation rate >0.002 nucleotides/generation, the proportion of bacteria expressing endonuclease below the threshold was dramatically increased and did not differ between selection and drift whereas in a population of 10000 bacteria, selection was the prime mover in generating bacteria below the endonuclease and arose after the G=40 threshold for established killswitches including the ‘Essentializer’.
One drawn conclusion from this modelling was the sensitivity of the current thanogen in the ccdA-ccdB binding and protease activity, likely on account of the higher likelihood of mutations occurring within the CDS. However, the failproof design of the killswitch ensures that the default state of the cell be cell death: dysfunction in protease, ccdB and endonuclease must arise before cell stasis/death occur. This, in part, justified the design of a more complex two-toxin system to supplant a single-plasmid T/A approach.
Ay, Ahmet, and David N. Arnosti. 2011. “Mathematical Modeling of Gene Expression: A Guide for the Perplexed Biologist.” Critical Reviews in Biochemistry and Molecular Biology 46 (2): 137–51. https://doi.org/10.3109/10409238.2011.556597.
Bethke, Qusheng Jin and Craig M. 2007. “The thermodynamics and kinetics of microbial metabolism.” American Journal of Science 643-677.
Chan, Clement T.Y., Jeong Wook Lee, D. Ewen Cameron, Caleb J. Bashor, and James J. Collins. 2016. “‘Deadman’ and ‘Passcode’ Microbial Kill Switches for Bacterial Containment.” Nature Chemical Biology 12 (2): 82–86. https://doi.org/10.1038/nchembio.1979.
Chandran, D., W. B. Copeland, S. C. Sleight, and H. M. Sauro. 2008. “Mathematical Modeling and Synthetic Biology.” Drug Discovery Today: Disease Models 5 (4): 299–309. https://doi.org/10.1016/j.ddmod.2009.07.002.
Cohen, Barner. 1954. “Studies on unbalanced growth in Escherichia coli.” Proc Natl Acad Sci U S A 885–893.
CONTOIS, D E. 1959. “Kinetics of bacterial growth: relationship between population density and specific growth rate of continuous cultures.” J Gen Microbiol 40-50.
D A Ratkowsky, R K Lowry, T A McMeekin, A N Stokes, and R E Chandler. 1983. “Model for bacterial culture growth rate throughout the entire biokinetic temperature range.” J Bacteriol 1222–1226.
D S Kompala, D Ramkrishna, N B Jansen, G T Tsao. 1986. “Investigation of bacterial growth on mixed substrates: experimental evaluation of cybernetic models.” Biotechnol Bioeng 1044-1055.
Fange, David, Harriet Mellenius, Patrick P. Dennis, and Måns Ehrenberg. 2014. “Thermodynamic Modeling of Variations in the Rate of RNA Chain Elongation of E. Coli Rrn Operons.” Biophysical Journal 106 (1): 55–64. https://doi.org/10.1016/j.bpj.2013.11.4487.
Frank, Heiner Hoffman and Michael E. 1965. “Synchrony of Division in Clonal Microcolonies of Escherichia coli.” J Bacteriol. 2 (89): 513–517.
Grant, Waclaw, Allen, Cicuta. 2014. “The role of mechanical forces in the planar to bulk transition in growing E. coli microcolonies.”
Hao, Nan, Adam C. Palmer, Alexandra Ahlgren-Berg, Keith E. Shearwin, and Ian B. Dodd. 2016. “The Role of Repressor Kinetics in Relief of Transcriptional Interference between Convergent Promoters.” Nucleic Acids Research 44 (14): 6625–38. https://doi.org/10.1093/nar/gkw600.
Hedges, Whitehead and. 2005. “Photodegradation and photosensitization of mycosporine-like amino acids.” Journal of Photochemistry and Photobiology B: Biology 115-121.
Jaime Sánchez-Navarrete, Nancy Jannet Ruiz-Pérez, Armando Guerra-Trejo & Julia Dolores Toscano-Garibay. 2020. “Simplified modeling of E. coli mortality after genome damage induced by UV-C light exposure.” Sci Rep 11240.
K A Presser, D A Ratkowsky, T Ross. 1997. “Modelling the growth rate of Escherichia coli as a function of pH and lactic acid concentration.” American Society for Microbiology Journals 2355-2360.
Leppek, Kathrin, Rhiju Das, and Maria Barna. 2018. “Functional 5′ UTR MRNA Structures in Eukaryotic Translation Regulation and How to Find Them.” Nature Reviews Molecular Cell Biology 19 (3): 158–74. https://doi.org/10.1038/nrm.2017.103.
Li, Mingji, Junshu Wang, Yanping Geng, Yikui Li, Qian Wang, Quanfeng Liang, and Qingsheng Qi. 2012. “A Strategy of Gene Overexpression Based on Tandem Repetitive Promoters in Escherichia Coli.” Microbial Cell Factories 11: 1–10. https://doi.org/10.1186/1475-2859-11-19.
Murray, James D. 2007. Mathematical Biology: I. An Introduction. Springer Science & Business Media.
n.a. 2013. “Reference Solar Spectral Irradiance: Air Mass 1.5.” National Renewable Energy Laboratory.
Vandervelde, Alexandra, Igor Drobnak, San Hadži, Yann G.J. Sterckx, Thomas Welte, Henri De Greve, Daniel Charlier, Rouslan Efremov, Remy Loris, and Jurij Lah. 2017. “Molecular Mechanism Governing Ratio-Dependent Transcription Regulation in the CcdAB Operon.” Nucleic Acids Research 45 (6): 2937–50. https://doi.org/10.1093/nar/gkx108.
Waclaw, Rosalind J Allen and Bartłomiej. 2018. “Bacterial growth: a statistical physicist’s guide.” Rep Prog Phys 016601.