Our Approach
We used two mathematical methods to design our microbial cell factories for increased perillyl alcohol(PA) yield. Our first model was a dynamic protein-oriented model to analyse the impact of different gene circuits with and without feedback control (for example testing the benefits of the efflux pump, engineered feedback system and tuning rates of different protein transcriptions).
Our second model was a linear genome-scale metabolic model of our E. coli . We used this framework to investigate gene editing options that could be used to increase our yield and consistency of our yield of PA prior to the addition of the novel efflux pumps or transcription feedback.
Our Dynamic Model
The Pathway Model
Our first model is an extension of the protein-based dynamic models developed
here. The core host model uses ODEs to model the transcription and translation of proteins used in the production of transport proteins, ribosomal enzymes, housekeeping proteins and enzymes dependent on the relative rates of action and availability of resources like energy (ATP) and free ribosomes.
We extended the model to include the MVA pathway in several stages: we first constructed a reduced version of our pathway maintaining linearity of reactions and relevant toxic metabolites along with the Michaelis-Menten style kinetic ODEs to represent the change in concentration of each of our metabolites depending on enzyme concentration. To find the Michaelis-Menten constants for our reactions we used data from here, the Brenda Database and approximations where the data wasn't available.
We then added the relevant enzymes to our model, this involved creating ODEs for the trancription and translation of the proteins used to make the proteins. To find appropriate rates for the activity needed to create these enzymes we found the amino acid length of the proteins involved using benchling or values from the original model.
Our next stage was to add in the toxicity of the relevant metabolites(s2, s3, s5, s6, s7) in our pathway, to do this we first created a Hill-type function-involving the concentrations of the toxic metabolites - that ranged from zero (modelling no toxic effect) to 1 (which would indicate maximal toxicity.) We then added a (1-toxic effect) multiplier to the conversion of our internal substrate to energy(ATP). We propose this model as the literature we reviewed proposed that the intermediate metabolites are toxic as they cause disruption (or perforation) of the E. coli cell membrane which hinders the action of ATP synthase. To make the effect of the toxicity more realistic we also introduced a death rate to the population growth so that the ability for the E. coli cells to live was dependent on the energy available to them, allowing our toxin to kill cells whereas the original model would only prevent cells from reproducing.
Given the complexity of our model we implemented the ODEs in MATLAB and used the matlab genetic algorithm tool to identify the rates of transcription for each of our enzymes that optimise total PA productionso we can investigate the effect of our efflux pumps and transcription rates assuming at improving the greatest rate of PA production.
Modelling Our Additions
To ensure we can test each of our constructs separately, we created a modular script that could run the pathway, the pathway with the efflux pumps, the pathway with the transcription feed back system and the pathway with both the efflux pump and the transcription feedback system so that the same script could be used to test and compare the various imposable conditions.
We first added Efflux pumps to our model, we did this by modelling the manufacturing of the efflux pumps to mirror the other transport proteins in our model, again we were able to find the length of the relevant proteins using benchling. In order to integrate our Efflux pumps into our model we modified our previous pathway to have both internal(inracellular) and external(extracellular) PA and limonene. This was an assumption of how the Efflux pumps would act as the Efflux pumps are specific to limonene but limonene and PA have similar enough chemical structures and interactions that the Efflux pump will transport them both across the cell membrane.
We used Michaelis-Menten style ODEs to convert the internal Limonene and PA to their external counterparts, for this reaction we modelled the reaction as having two competing substrateswith the efflux pumps being modelled as the enzyme catalysing these reactions.
Our tests comparing a batch culture of E. coli with and a batch culture without the Efflux pumps showed very promising improvements for both total PA yield and limonene yield for the cells expressing Efflux pumps to the cells that don't, making it clear that the efflux pumps were worth pursuing.
The other option we added to our model was the transcription feedback loop, we trialled modelling this change to the rate of transcription of the first enzyme in our model (in reality this enzyme represents 3 enzymes in the actual E-coli MVA pathway) using several different techniques: the first approach was to have the transcription rate directly altered by the toxic effect as a (1-alpha*toxm) multiplier where alpha ranges between 0 and 1 but as this didn't include the envelope pathway response delay it caused changes that were unrealistically rapid; we therefore changed our approach and added the Tet repressor (TetR) proteins to our model to act as intermediataries between the toxic effect and transcription rate change - a more realistic,mechanistic approach.
To model the affect of the TetR protein concentration on the transcription rate of our first enzyme (E1), we used a Hill function to model the binding of the proteins to the relevant promoter regions that control the transcription rate, we refer to the value of this Hill function as the Tet factor. We then added a (1- Tet factor) multiplier to the transcription rate of the relevant enzyme and ran tests to compare the performance of batches with and without the addition of our feedback loop on our E-coli without the Efflux pumps present.
Even with the Tet factor Hill function weighted closer to 1, the transcription feedback system had a significantly weaker effect on the amount of PA being produced than the efflux pump and although the feedback system slightly increases the yield of limonene, it decreases the yield of PA. We then added artificial degradation to the E1 enzyme so that the concentration of E1 would be more responsive to the TetR.
This responsiveness can be easily seen in the plateau stage of the graph above. It was clear from our model that the transcription feedback seemed to decrease rather than increase our PA yield. Tis shows the system needs significant tuning; a process we leave for future work.
Our Metabolic Model
Secondly we took a `conventional' metabolic engineering approach and utilsied a genome-scale model of E. coli to identify gene knockouts which would improve the yields of PA. We extended the genome-scale model to include efflux pumps but not the dynamic non-linear transcription feedback system.
Our goal with this model was to utilise a series of algorithms to run OptForce, a MATLAB Cobra Toolbox function that is used recommend upregulations,downregulations or genetic deletions to a subset of the genes in the inputted genome. To do this we used the iB21_1397 BIGG model: this model consists of a stoichiometry matrix of coefficients of the 2741 typical reactions for the BL21 E-coli strain our team is genetically engineering.
We added our extended MVA pathway reactions to the model, adding the efflux pump as two export reactions of limonene and PA out of the cell. For this purpose, we used the full detailed metabolic pathway shown below:
This linear model can be visualised as a flow network of the carbon atoms of glucose through the cell with extracellular glucose as a source node and and extracellular PA and biomass as the sink nodes of interest. Our first algorithm was Flux Balance Analysis (FBA), a form of cut-set analysis that we used to separately find the maximal output of PA and biomass from our network along with the vectors of fluxes that allow for the maximum flux of carbon through the biomass reaction or the PA export reaction.
We then compared the multiple optimal vectors in an algorithm called Flux Variance Analysis (FVA), which utilises the degeneracy of the optimal flux solution vectors. It compares the overlap of the range of the flux vectors that maximise biomass and the flux vectors that maximise the export of PA to create MUST sets to suggest up and down regulations for each reaction and the related genes.
These MUST sets were then used in the OptForce algorithm, which computed the minimal set (FORCE set) of genetic modifications that would be needed to meet the desired overproduction yield of PA.
Unfortunately OptForce was unable to find these regulations in a reasonable timeframe. Therefore we tried a new tool gcfront . This tool seeks to identify metabolic engineering strategies where the production of the metabolite of interest is coupled to biomass production processes.
gcFront outputted 37 pareto optimal designs but even with 4 knockouts the maximum PA flux was only 0.01272 when knocking out the B21_03847, B21_00676 and B21_02636 genes, a small proportion of the 2.3070 mmol/hour that FBA says the E-coli should be capable of. Even more poorly, the optKnock algorithm outputted three reactions the deletion of which affected neither the biomass or PA production when deleted.
These results suggest there is a hard trade-off between growth and PA production which causes there to be few, or potentially no good gene knockout strategies that allow PA synthesis whilst also creating a significant amount of growth.