Team:William and Mary/Model
Modeling
.
Overview
We constructed and implemented a complex, novel mathematical model to answer the following three questions:
- Does our model indicate feasibility of TheraPUFA in decreasing viral load using parameters sourced from the literature?
- Is there a way to optimize parameters such that the model indicates feasibility?
- If we incorporate diffusion in a spatial environment and transcriptional stochasticity into our model, is TheraPUFA still feasible?
Using the specific set of parameters that we sourced from literature, our model initially indicated that TheraPUFA was not feasible as an antiviral drug. However, we utilized our model to optimize these parameters and identify values under which TheraPUFA was effective in significantly decreasing viral load. With the optimized parameters, TheraPUFA was still effective even when transcriptional stochasticity and spatial aspects were accounted for.
.
Modeling PUFA Synthesis & Export
To model the probiotic’s effect on viral load, we first modeled the probiotic’s synthesis and export of AA and DHA. The AA and DHA “synthesis and export” models constitute two components of a larger model that simulates the interaction of PUFA with viruses and the natural clearing of PUFA from the nasal cavity due to mucociliary clearance.
As described in detail on our Design page, we designed multiple systems for PUFA export in order to explore the possibility of implementing TheraPUFA in a gram positive chassis. The finalized model utilizes the systems designed to export AA and DHA from our gram negative chassis N. cinerea, since our cytokine sensor contains both outer and inner membrane components and thus requires expression in a gram negative chassis. However, we still designed two systems to export AA and one system to export DHA from a gram positive chassis. A graphic summary of all the export systems we designed is shown below. “+” and “-” refer to expression in a gram positive or negative chassis, respectively. The production and export of PUFA under each system was modeled according to Michaelis-Menten kinetics.
Michaelis-Menten Rate Law
In our circuit model, the majority of steps for the production and export of fatty acids were modeled through the use of the Michaelis-Menten rate law, which used a quasi-steady-state assumption to simplify the following enzymatic reaction process:
This was converted into a quicker and more usable form which described the reaction rate of a substrate becoming a product in the following way:
Where Vmax is equal to the maximum reaction rate based on the amount of enzyme present and Km is the amount of substrate needed in order to achieve half of the maximum reaction rate. In our model, Km and Vmax were both normalized into the same units through dimensional analysis, such that the equations simulated nanomoles of PUFA per single cell.
Vmax can be converted to units of nmol·min-1·cell-1. To accomplish this conversion, the metric of copies of enzyme per typical cell was estimated from the literature. Next, the mass of each enzyme in kilodaltons was found and converted to mg. By multiplying the copies of enzyme in the cell by the mass of each enzyme in mg, the total mg of enzyme per cell can be obtained. This mass can be multiplied by Vmax in nmol·min-1·mg-1 to yield Vmax in nmol·min-1·cell-1.
As a concentration, Km is expressed in units of uM (micromolar, or micromoles per liter volume). Since all other metrics in the model are expressed in nanomoles per cell rather than micromoles per liter, Km was also converted to nanomoles per cell. An E. coli cell is about 1 um3, and 1 mL or 1 cm3 = 1012 um3.
In the parameter spreadsheet provided under "PUFA Synthesis and Export Models," parameter tables are provided which state the Vmax and Km values as they were found in the literature, as well as adjusted values (intermediate conversions not shown).
Model Assumptions
- Starting substrates (acetyl-CoA and malonyl-CoA) are in abundance and do not change over time. In gram negative cells, we will assume that Coenzyme A, which "activates" fatty acids into fatty acyl CoA, is also in abundance.
- Fatty acids are not degraded after synthesis. (This would involve FadE knockout in gram negative chassis. In some species that lack β-oxidation pathways, such as gram positive S. aureus, no additional engineering is necessary to prevent degradation).
- Coenzyme A is in abundance and does not limit activation.
- Once released into the periplasm, fatty acids do not move “backwards” into the cytoplasm. We can assume that we knock out FadD, an inner-membrane associated protein that uses ATP to move protonated free fatty acids from the periplasm into the cytoplasm.
- Once released into the extracellular environment, DHA/AA does not re-enter the cell through the outer membrane. For example, in gram negative species, we can assume we have knocked out FadL, which transports fatty acids into the periplasm from the extracellular environment.
- We will assume that the concentration of DHA/AA outside of the cell is always 0 due to diffusion of PUFA and robust mucociliary clearance. Because of this assumption, we can also assume that free DHA/AA outside the cell will not inhibit the growth of the cell (free fatty acids can inhibit bacterial growth, but we are assuming that mucociliary clearance, etc., will carry PUFA away from cells quickly).
PUFA Synthesis and Export Models
To view the ordinary differential equations for each designed PUFA export system, simply click on the export system name. Only the first two export models listed (AA and DHA export in a gram negative chassis) were implemented in the overall model to simulate TheraPUFA's effect on viral load. The remaining three models listed are alternatives that we designed to explore implementing TheraPUFA in a gram positive chassis.
A spreadsheet containing parameters utilized in the models described below can be found here.
- AA is synthesized by PUFA synthase and incorporated into phospholipids of the inner cell membrane. Unlike DHA, which appears to initially accumulate as intracellular free fatty acid, AA appears to accumulate within phospholipids—no additional engineering is required to divert synthesized AA to phospholipid synthesis (Metz et al., 2009, Ujihara et al., 2014). Therefore, when modeling AA production, we collapse synthesis and phospholipid incorporation into one step.
- cPLA2 localized to the periplasm utilizing a CycA' localization tag from Neisseria spp. releases free AA from the inner membrane (Turner et al., 2005; Rosa & Rapoport, 2009).
- Exit from the periplasm occurs through the outer membrane through a Tol-C independent mechanism (Tong et al., 2019). It is unlikely for fatty acids to passively diffuse through the outer membrane, especially considering the impermeability of bacterial membranes to long chain fatty acids (Schwenk et al., 2010; Kamp & Hamilton, 2006). So, we will assume a protein-mediated diffusion mechanism rather than simple diffusion. Specifically, we will assume that a FadL-like protein allows periplasmic free AA through the outer membrane and into the extracellular environment. FadL is a long-chain fatty acid transporter that moves free fatty acid into the periplasm (van den Berg et al., 2004). We will assume that there is a similar protein that moves free fatty acid in the opposite direction.
These steps can be represented in the following equations:
α(t) represents the generation of AA within membrane phospholipids, or step 1. As mentioned previously, the production of AA and its incorporation in phospholipids is condensed in α(t).
β(x1) represents cleavage of AA from membrane phospholipids by cPLA2, a phospholipase with a preference for AA, and release into the periplasm, or step 2.
γ(x2) represents exit of free AA from the periplasm and into the extracellular environment, or step 3.
Accounting for the parameters in the spreadsheet provided, the equations become:
All measurements above have been scaled by 1013 and thus represent nmol x 1013 per cell per minute.
- DHA is synthesized as a free fatty acid by a Schizochytrium PUFA synthase (Metz et al., 2009). An additional phosphopantetheinyl transferase from Nostoc cyanobacteria is necessary for DHA synthesis when the DHA synthase is expressed heterologously.
- Free DHA is activated into DHA-CoA by acyl-CoA synthetase lacsA from diatom Thalassiosira pseudonana (Tonon et al., 2005).
- Activated DHA is incorporated into membrane phospholipids and inserted into the membrane by enzymes (PlsC incorporates fatty acids into the sn-2 position. We will only consider PlsC for this step, because potential chassis Neisseria cinerea utilizes the PlsX/Y/C system and lacks a PlsB enzyme capable of using DHA-CoA) (Sohlenkamp & Geiger, 2016).
- DHA is released from the inner membrane and into the periplasm by a phospholipase. Assume iPLA2, a human phospholipase with a preference for DHA, has been localized to the periplasm using a CycA’ localization tag endogenous to Neisseria spp. (Rosa & Rapoport, 2009; Turner et al., 2005).
- Exit from the periplasm occurs through the outer membrane through a Tol-C independent mechanism (Tong et al., 2019). It is unlikely for fatty acids to passively diffuse through the outer membrane, especially considering the impermeability of bacterial membranes to long chain fatty acids (Schwenk et al., 2010; Kamp & Hamilton, 2006). So, we will assume a protein-mediated diffusion mechanism rather than simple diffusion. Specifically, we will assume that a FadL-like protein allows periplasmic free DHA through the outer membrane and into the extracellular environment. FadL is a long-chain fatty acid transporter that moves free fatty acids into the periplasm (van den Berg et al., 2004). We will assume that there is a similar protein that moves free fatty acids in the opposite direction.
The steps described above can be modeled with the following equations:
α(t) represents generation of DHA as an intracellular free fatty acid, or step 1.
β(x1) represents activation of DHA to DHA-CoA by the acyl-CoA synthetase.
γ(x2) represents incorporation of activated DHA into membrane phospholipids by the endogenous enzyme PlsC.
μ(x3) represents cleavage of DHA from phospholipids and release into the periplasm by iPLA2, a phospholipase with specificity for DHA.
κ(x4) represents the diffusion of DHA from the periplasm as mediated by a FadL-like protein.
Accounting for the parameters described above, the final equations used to model the synthesis and export of DHA from a gram negative chassis are:
All measurements above have been scaled by 1013 and thus represent nmol x 1013 per cell per minute.
- AA is synthesized by the A. marina PUFA synthase (Ujihara et al., 2014).
- AA diffuses to the inner leaflet of the cell membrane, where it embeds itself. The FakA/B3 complex phosphorylates the embedded AA to AA-PO4 (Gullet et al., 2019).
- AA-PO4 is utilized by PlsY and inserted into the cell membrane as the sn-1 component of phospholipids (Lu et al., 2006).
- Thraustochytrid 145138, a phospholipase capable of cleaving at the sn-1 position, cleaves phospholipids and releases free AA into the gram positive “periplasm” (Ishibashi et al., 2019).
- Free AA diffuses freely/passively through the cell wall.
The steps described above can be modeled with the following equations:
α(t) represents generation of AA as an intracellular free fatty acid, or step 1.
β(x1) represents phosphorylation of AA to AA-PO4 by FakB3.
γ(x2) represents incorporation of phosphorylated AA into membrane phospholipids by PlsY.
μ(x3) represents cleavage of AA from membrane phospholipids by Thraustochytrid phospholipase 145138.
Accounting for parameters, the final equations used to model the synthesis and phospholipase-mediated export of AA from a gram positive chassis are:
All measurements above have been scaled by 1013 and thus represent nmol x 1013 per cell per minute.
Alternatively, it is possible that the PlsC enzyme endogenous to our potential gram-positive chassis can act upon AA attached to the acyl carrier protein (ACP) of the PUFA synthase. In this case, incorporation of AA into membrane phospholipids could occur without additional engineering, and AA would be inserted at the sn-2 position of phospholipids, necessitating a PLA2 phospholipase rather than a PLA1.
Whereas DHA accumulated extracellularly in E. coli as a free fatty acid, AA was incorporated into membrane phospholipids (Metz et al., 2009). Without additional research on heterologous expression of AA synthases within gram positive bacteria, it is unclear whether gram positive cells can incorporate the PUFA within their phospholipids.
This process could be modeled with the following equations, where x1 represents the change in AA bound to ACP of the AA synthase, x2 represents the change in phospholipid-incorporated AA, and x3 represents the change in free, secreted AA in the extracellular environment:
α(t) represents generation of AA by the AA synthase.
β(x1) represents the removal of AA from ACP and subsequent insertion into membrane phospholipids by PlsC.
γ(x2) represents the release of AA from membrane phospholipids by cPLA2.
Accounting for the parameters in the spreadsheet provided, the equations become:
All measurements above have been scaled by 1013 and thus represent nmol x 1013 per cell per minute.
- DHA is synthesized as a free fatty acid by a Schizochytrium PUFA synthase (Metz et al., 2009). An additional phosphopantetheinyl transferase from Nostoc cyanobacteria is necessary for DHA synthesis when the DHA synthase is expressed heterologously.
- DHA diffuses to the inner leaflet of the cell membrane, where it embeds itself. The FakA/B3 complex phosphorylates the embedded DHA to DHA-PO4 (Gullet et al., 2019).
- DHA-PO4 is utilized by PlsY and inserted into the cell membrane as the sn-1 component of phospholipids (Lu et al., 2006).
- Thraustochytrid 145138, a phospholipase capable of cleaving at the sn-1 position, cleaves phospholipids and releases free DHA into the gram positive “periplasm” (Ishibashi et al., 2019).
- Free DHA diffuses freely/passively through the cell wall.
These steps can be represented in the following equations:
α(t) represents generation of free, intracellular DHA by the DHA synthase.
β(x1) represents the phosphorylation of DHA into DHA-PO4 by FakB3.
γ(x2) represents incorporation of phosphorylated DHA into membrane phospholipids by PlsY.
μ(x3) represents cleavage of DHA from phospholipids by the Thraustochytrid phospholipase 145138.
With parameters accounted for, the final equations are:
All measurements above have been scaled by 1013 and thus represent nmol x 1013 per cell per minute.
In the absence of a β-oxidation pathway to degrade fatty acids, pathogenic bacteria such as S. aureus have two options to detoxify PUFAs secreted by the human body: 1) incorporate PUFAs into their cell membranes, or 2) efflux PUFAs. To efflux linoleic acid (LA) and AA, S. aureus utilizes an efflux pump called FarE, regulated by FarR (Alnaseri et al., 2015; Alnaseri et al., 2019). By removing FarE from the control of FarR and placing it beneath the control of an inducible promoter (such as the pspA promoter), we can make FarE responsive to cytokine concentrations.
- AA is synthesized utilizing an A. marina PUFA synthase. E. coli seems capable of activating free AA and inserting it into the membrane without further engineering (Ujihara et al., 2014). While it is unclear whether gram positive bacteria would similarly insert AA into membranes, PlsX may be able to utilize AA-ACP as a substrate, converting it to AA-PO4, which is then used for phospholipid synthesis by PlsY. PlsC may also be able to utilize AA-ACP as a substrate (Lu et al., 2006). In our model, we assume our gram positive chassis can utilize AA for membrane synthesis, and we account for a decrease in intracellular free AA that can be effluxed by FarE. In other words, the amount of free AA available is equal to AA synthesized by the PUFA synthase minus AA used by PlsX and AA used by PlsC.
- Free AA is pumped from the cell by rhe FarE efflux pump (Alnaseri et al., 2015; Alnaseri et al., 2019).
These steps can be represented in the following equations:
α(t) represents generation of AA attached to ACP on the AA synthase.
β(x1) represents the utilization of AA by PlsX.
With parameters input, the equations become:
All measurements above have been scaled by 1013 and thus represent nmol x 1013 per cell per minute.
.
Modeling TheraPUFA's Antiviral Effect in Human Nasal Cavity
To model the effect of TheraPUFA on viral load, we altered the target cell limited model (Zitzmann & Kaderali, 2018) to fit the nasal cavity, incorporated the effect of DHA and AA on the viral titer, and added the production of cytokines.
Variables
I: infected cells
V: viral load
T: concentration of TNF-α
F: concentration of IFN-γ
Parameters
d : death rate of susceptible cells
k : viral infection rate
δ : death rate of infected cells
p : rate for infected cells to produce virus
c : viral clearance rate by all mechanisms
εpAD : efficiency for AA and DHA work together to prohibit the viral replication
εTA : the efficiency of AA to promote the production of TNF-α
εTD : the efficiency of DHA to prohibit the production of TNF-α
εFAD : the efficiency of AA and DHA work together to prohibit the production of IFN-γ
pT : the production rate of TNF-α by infected cells
pF : the production rate of IFN-γ by infected cells
dT : decay rate of TNF-α
dF : decay rate of IFN-γ
γ : strength of immunosuppression
Click here to download our model parameters
Note: How to get all ε terms?
Based on literature (Best, 2018; Goyal, 2020), the antiviral effect of a single drug is usually defined to be:
where C indicates the concentration of drug, EC50 indicates the concentration required to achieve half of the maximum effect.
Based on the literature (Kouizumi,2014), the efficiency of combination dose of drug A
and B,if they are additive, are expressed as the following equation:
To adapt this into our model, for all epsilon terms where only one substance adds effect, we subtituate C with our own variable into equation 1. When two substances work together to add an effect, we substitute C*A and C*B correspondingly.
Epsilon terms are bridges between the in the nose target cell limited model and PUFA production and export model. Since the efficiency is calculated based on the concentration of AA and DHA, which is dependent on the result from the PUFA synthesis and export model.
Equations
Explanation
Equation 1 accounts for the rate of change of susceptible cells (S). Susceptible cells are produced at a constant rate λ, and they die at a constant rate, d. When a susceptible cell gets infected, it becomes an infected cell (I) and is not considered a susceptible cell anymore. In the first equation, we can see that the amount of susceptible cells increases through production, and decreases through cell apoptosis or infection.
Equation 2 accounts for the rate of change of infected cells. The amount increases when susceptible cells are converted to infected cells through infection. Hence, the increasing term should be the same as the one accounting for infection process in equation one but with opposite sign. Additionally, infected cells die at a constant rate δ.
Equation 3 accounts for the rate of change of viral load. It gets increased when infected cells produce viral particles at a constant rate, which normally should be p. However, since AA and DHA prohibit the production process, this rate is reduced to (1 - εpAD). The viral load decreases through clearance by all mechanisms (innate immune response, adaptive immune response, etc.). Here we are assuming optimistically that we are able to use a single parameter c to cover all of those.
Equation 4 accounts for the rate of change of TNF-α. Infected cells produce TNF-α at rate pc, which can probably be promoted by AA and prohibited by DHA. Cytokines decay at rate dc.
Equation 5 is exactly the same with equation 4 except that AA and DHA both have a prohibitive effect on the production of IFN-γ.
.
Incorporating Spatial Aspects
The most common approach in modeling viral infections and the effect of antiviral treatments is through a system of differential equations, we have decided to go beyond this in incorporating the spatial aspects of the nasal cavity by implementing a 10X10 grid in which each cell of our grid contains its own set of our system of differential equation with its own initial values. Thus the concentration of each value is tracked independently for each cell. This accurately represents how bacterial cells naturally cluster in the nasal cavity. To account for non-homogeneous distribution, different initial conditions are used for each cell. For example, if we want a gradient style initial distribution for viral particles, we assign more viral particles in the cells to the left compared to those to the right initially.
Dividing space into grids can also account for the diffusion of small particles, namely AA, DHA, TNF-α, IFN-γ, and virus. The mathematical way to express that is through adding a diffusion term to equations for those variables in each cell. At every time step, the concentration of each grid cell is compared with the sum of the neighboring 8 grid cells, that is unless the cell exists on a boundary position of the grid. The difference in concentration between the cells determines the degree to which the particles diffuse. This difference is then multiplied by the diffusion constant (expressed as θ in the equations below).
By incorporating spatial aspects in our model, we have gone above and beyond the normal Target cell limited model in modeling viral infections.
.
Testing Feasibility
- Our first main objective of creating and running this model is to determine the feasibility of the design of our probiotic. Through modeling, we were able to quantify the peak viral load that exists without the production of polyunsaturated fatty acids (PUFAs) arachidonic acid (AA) and docosahexaenoic acid (DHA). We were then able to test if the synthesis and export of the PUFAs by our probiotic can cause a significant decrease in viral load by orders of magnitude.
- Our second main objective is to quantify the improvements we would need to make our design feasible.
First, when we ran the model without introducing synthesis and export of PUFAs, and we found the peak viral load to be 1.71e10 RNA copies/mL. Then we ran the model after incorporating PUFA synthesis and export using parameters we found in the literature and used an initial probiotic cell value of 1e8 without incorporating any replenishment of the probiotic. This resulted in a decrease of only 8.24e4 RNA copies/mL from the peak viral load, which is not significant at all.
However, it was too soon to conclude that our design was not feasible. As there was some ambiguity and variety in the parameters sourced from the literature, thus the implementation of PUFA synthesis and export did not have to be limited to the values initially used. Its effect on viral load could be increased through various means, such as the replenishment of the probiotic at reasonable intervals and the clearance rate of probiotic and PUFA from the nasal cavity. Below is a list of parameters that we have varied in order to examine their effect on peak viral load. Our strategy was to keep all other conditions the same, and to vary only one parameter at a time (using several values spanning a wide range). Note that since we aimed to find patterns in the relationships between parameters and difference from peak viral load, we have used some values that may seem to be biologically inaccurate.
| # | ||||
|---|---|---|---|---|
| 1 | Number of probiotic cells administered | Increase | 1e8, 1e9, 1e10, 1e11 | This value can be increased by increasing the value of the each dose |
| 2 | Interval for probiotic replenishment (probiotic is administered once every i hours) | Decrease | No replenish, Every 1, 2, 4, 6, 8, 10, 12, 16, above 16 hours | By replenishing more frequently |
| 3 | Probiotic clearance rate | Decrease | (in terms of fold of the original value) 1, 5e-1, 1e-1, 5e-2, 1e-2 | By adding adhesive gel / polymer matrix to increase attachment of probiotic to epithelial surface. |
| 4 | TNF-α sensor sensitivity (concentration of TNF-α above which DHA start to be produced) | Decrease | Above 20, 20, 18, 16, 14, 12, 10, 8, 6, 4, 2, 1 | By choosing sensor that can trigger the switch to production from AA to DHA at a lower TNF-α concentration, and / or by improving our existing sensor with synbio techniques |
| 5 | IFN-γ sensor sensitivity (concentration of IFN-γ above which DHA start to be produced) | Decrease | Above 2.6, 2.6, 2.5, 2.4, 2.2, 2.0, 1.8, 1.6, 1.4, 1.2, 1.0 | By choosing sensor that can trigger the switch to production from AA to DHA at a lower IFN-γ concentration, and / or by improving our existing sensor with synbio techniques |
| 6 | Maximum value for TNF-α | Increase | 0,100,150,200,220,300,400,500,600,700 | see if our circuit helps during a severe cytokine storm |
| 7 | Maximum value for IFN-γ | Increase | 0,100,200,300,400,500,600,700 | see if our circuit helps during a severe cytokine storm |
| 8 | Gap time (time required for DHA/AA to start having an effect) | Decrease | Every combination of the replenish time equals 2, 4, 8, 16, 24 and gap time from 0 to 48 with interval equals 2. | In the literature the effect on viral load is started to be measured 24 hours after the implementation, so the gap time can actually be anytime before that. |
| 9 | Production rate of AA/DHA | Increase | (in terms of fold of the original value) 1, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50 | Through metabolic engineering and synbio optimizations. For instance, we can delete fatty acid synthesis genes native to our chassis to increase PUFA |
| 10 | Decay rate of AA/DHA | Decrease | (in terms of fold of the original value) 1, 0.75, 0.7, 0.6, 0.5, 0.4, 0.25, 0.1, 0.07, 0.05, 0.03, 0.01 | Maybe there are ways we can increase PUFA incorporation rate into epithelial cells with some chemical additive or biological product to the nasal spray. |
Varying Initial Number of Probiotic Cells and Probiotic Replenishment Time Interval
Relationship Between Initial Number of Probiotic Cells and Difference in Peak Viral Load
After analyzing the three graphs shown above and the data points vertically, it is clear that when the initial number of probiotic cells has been increased by 10 fold, viral load difference is also increased by roughly 10 times for the same replenishment time interval. This relationship is also reflected by the 10-fold increase for the fitted curve when the initial number of probiotic cells is increased by 10 times. Therefore, there seems to be a linear relationship between the initial number of probiotic cells and the difference in peak viral load (when PUFA production is not introduced compared to when it is). When relating this relationship back to the administration of our probiotic, we should aim to choose the highest number of probiotic cells per that is possible under biological and safety constraints.
Logarithmic Relationship Between the Difference in Peak Viral Load and PUFA Production
The three graphs above show that the replenishment time interval (number of hours between replenishment of the probiotic to return the number of probiotic cells back to the initial value) has a logarithmic relationship with the difference in peak viral load (when PUFA production is not introduced compared to when it is, as shown in the graphs above). Another observation we have made is that when we replenish the probiotic less than one time every 16 hours, there is no difference in peak viral load from when there is no replenishment incorporated at all. When relating this relationship back to the administration of our probiotic, we should aim to replenish the probiotic as frequently as possible under biological and safety constraints.
Varying Probiotic Clearance Rate
Varying Sensor Sensitivity to TNF-α
Observations Regarding Sensor Sensitivity to TNF-α
Given that TNF-α concentration only reaches a maximum concentration of around 20 pM, when the sensor sensitivity is above 20 pM the concentration level needed to activate the sensor and trigger the production of DHA is never achieved. Therefore, the results are no different from the results we found when sensitivity to IFN-γ above this threshold. When sensor sensitivity to TNF-α becomes less than or equal to 20 pM, we begin to see a decrease in peak viral concentration.
While the version where AA is consistently produced outperforms the version where AA and DHA are switched, the improvement is not very significant. Thus it is important to consider the beneficial effects which switching to DHA has for the inflammatory response inside the body.
We note that change in sensitivity is not vastly different from the value sourced in the literature and that this is achievable through directed evolution of our sensor.
Varying Sensor Sensitivity to IFN-γ
Observations Regarding Sensor Sensitivity to IFN-γ
Given that IFN-γ concentration only reaches a maximum concentration of around 2.4 pM, when the sensor sensitivity is above 2.6 pM, the concentration level needed to activate the sensor and trigger the production of DHA is never achieved. Therefore, the results are no different from the results we found when sensitivity to IFN-γ above this threshold. When sensor sensitivity to IFN-γ becomes less than or equal to 2.6 pM, we begin to see a decrease in peak viral concentration.
In order to begin decreasing peak viral load, the sensor sensitivity needs to be as low as 2.6 pM. This is a significant change from its original where is 200 pM. These results show that IFN-γ may not be feasible as a trigger for our sensor. As this sensor can detect concentrations of either TNF-α or IFN-γ, we may aim to adjust the sensitivity for TNF-α alone.
Varying Maximum Concentration of TNF-α and IFN-γ
Observation Regarding Varying Maximum Cytokine Concentrations
When TNF-α’s maximum concentration is increased to be over 180 fold times more than that of the original maximum concentration, DHA production is activated, but the amount is too small to be shown on the graph. After TNF’s maximum concentration reaches above 200 fold times more than the original maximum concentration, the difference in peak viral load increases quickly, indicating a more significant decrease in peak viral load when compared to peak viral load without the incorporation of PUFA production. After increasing the maximum concentration of TNF-α by over 300 fold, the difference in peak viral load increases in a linear way. This is a similar trend seen with the varying of IFN-γ’s maximum concentration. However, the difference in peak viral load caused by changing TNF-α’s maximum concentration is greater than the one caused by changing IFN-γ’s maximum concentration.
If the sensitivity of the sensor is unchanged, It is clear that we would need to have a maximum concentration of TNF-α and IFN-γ higher than the average concentrations found in the nasal cavities of patients with respiratory viral infections in order to make a significant difference in peak viral load. It is important to note that these values were sourced from a nasal lavage and each paper quotes a different maximum concentration. In addition, it is entirely possible that the concentration of TNF-α and IFN-γ might in fact be higher in the nasal cavity.
Varying Maximum Concentration of Either TNF-α or IFN-γ
Observation Regarding Varying Either TNF-α or IFN-γ
The circuit design in which PUFA production switches between AA and DHA results in a smaller decrease in viral load due to sensor activation. Yet it is important to note that this result is not very significant. Thus it is important to consider the effect which switching from AA and DHA would have on the inflammatory response of nasal epithelial cells.
Varying Gap Time (Time Required for DHA/AA to Have an Effect)
Observation Regarding Varying Gap Time
Replenishment time of 16 hours and 24 hours result in similar outcomes, thus we have collapsed them in the graph above. This is due to the fact that peak viral load is achieved at 16 hours (i.e. the threshold time), once replenishment time is greater or equal to this time, changing replenishment time will have no influence. If replenishment time is smaller than the threshold time, then as replenishment time gets smaller, the viral load difference becomes greater, and peak viral load becomes lower. If gap time is greater than the threshold time, then there will be no effect on peak viral load, which results in a null effect on the viral load, as shown by the graph. If gap time is smaller than this threshold time, then as gap time becomes smaller, the viral load difference becomes greater. It can be explained by the following: a smaller gap time will cause the antiviral effects of AA/DHA to occur earlier, causing a greater difference in viral load as AA and DHA will have a longer time period for their antiviral effects to take place.
Varying Production Rate of AA and DHA
Observation Regarding Varying AA and DHA Production Rates
An increase in AA and DHA production results in an increase in peak viral load difference which follows a logarithmic relationship. Thus, the increase of AA and DHA would not yield a significant effect after a certain increase in fold (based on the graph it is roughly 20 fold).
Varying AA/DHA Clearance Rate
After increasing the production rate the clearance rate of AA and DHA is varied and its results are noted below.
Observation Regarding Varying AA/DHA Clearance Rate
As the PUFA decay rate becomes lower, the difference in viral load increases at a greater rate.
.
Results with Final Adjusted Parameters
After analyzing the effect of varying single parameters on the difference in peak viral load between the scenario in which PUFA production was not introduced into the model (1.71e10 RNA copies/mL) and in the scenario provided in figure 12, we vary them in combinations based on our observations. We define a significant effect as being 2-3 orders of magnitudes lower than the peak viral load in the scenario in which no PUFA production is incorporated. After testing 8 combinations of parameters, we found that the following set of parameters gives us this significant result. Specifically, peak viral load has decreased from 1.71e10 to 8.72e8. Additionally, the viral load also decreases at a faster rate and eventually reaches the value of 1e6. On the other hand, the viral load remains as high as 1e9 RNA copies/mL if there is no PUFA incorporated. This suggests successful antiviral effects of PUFA given this set of parameter adjustments.
Final Parameters
| Probiotic number | as great as possible in practical range | 2e11 |
| Interval for probiotic replenishment | as great as possible in practical range | Every 6 hrs |
| Probiotic clearance rate | a reasonable value after which the effect is not that significant (roughly 1e-1 based on the figure). | 1e-1 of the original value |
| TNF-α sensor sensitivity | Switch production of DHA and AA above threshold; sensor as sensitive as we can possibly get | 10 pM |
| IFN-γ sensor sensitivity | Leave unchanged since adjusting the TNF-α sensor alone helps to switch the DHA production on. | Original value |
| Gap time | Get greatest difference when it is 0 | 0 |
| Production rate of AA/DHA | More significant impact before 20 folds | 20 |
.
Robust to Higher Initial Viral Load?
Observation
As the initial viral load increases, the viral load difference increases dramatically. After reaching the peak the viral load difference begins to decrease linearly. Under this set of parameters, the peak is around 190 fold.
It is best for us to find the peak viral load and use that value for the initial viral load to lower the viral load peak. Based on this graph, we should possibly test cases where initial viral load goes to 400 and above. That’s where the difference starts to decrease, so we need to test if under this scenario our set of parameter choices is still able to have a relatively significant effect on viral load. However, for the sake of making our tests more robust, we have simulated scenarios when initial viral load has been increased to 200, 300, 400, 500 viral particles respectively.
Test 1: When initial viral load has been increased from 100 viral particles to 200 viral particles.
Peak viral load w/o PUFA: 1.72e10
Difference:1.61e10
Test 2: When initial viral load has been increased from 100 viral particles to 300 viral particles.
Peak viral load without PUFA production incorporated: 1.72e10
Difference in peak viral load:1.61e10
Test 3: When initial viral load has been increased from 100 viral particles to 400 viral particles.
Peak viral load without PUFA production incorporated: 1.73e10
Difference in peak viral load: 1.61e10
Test 4: When initial viral load has been increased from 100 viral particles to 500 viral particles.
Peak viral load without PUFA production incorporated: 1.73e10
Difference in Peak Viral Load: 1.6e10
.
Robust When Stochasticity is Incorporated?
Since the production and export of AA and DHA require several enzymes, it is necessary to also consider the fact that there is stochasticity involved in this process. In order to account for that, we add Gaussian white noise terms, with mean 0 and standard deviation 1 order of magnitude less than the original parameter, to all parameters involved in the production and export of AA and DHA equations. Then we test if our set of parameters is still robust.
Since the stochastic process leads to different results from trail to trail, we’ve tested several times and look for the average.
| w/o stochasticity | 1.71e10 | 8.72e8 | 1.62e10 |
| 1 | 8.77e8 | 1.62e10 | |
| 2 | 8.72e8 | 1.62e10 | |
| 3 | 8.74e8 | 1.62e10 | |
| 4 | 8.68e8 | 1.62e10 | |
| 5 | 8.66e8 | 1.62e10 | |
| Mean of the 5 trials | 8.7e8 | 1.62e10 |
As we can see from the bottom row, with stochasticity incorporated, our set of parameters is still able to decrease the peak viral load by the same order of magnitudes. Actually, the result is pretty approximate to the case where stochasticity is not incorporated at all. That’s because all parameters involved in the production and export of DHA and AA are on really small order of magnitudes. If we choose to add Gaussian white noise with standard deviation 1 order of magnitude smaller, there will not be too much fluctuations added.
.
Robust When Spatial Factors are Incorporated?
So far we have been simulated using only the ODE model. While it is the most commonly used strategy for most viral dynamics models, it assumes homogeneous distribution for all variables, which is usually not the case in reality. Specifically in our model, viral infection usually starts from one end of the nasal cavity. Therefore, what’s more approximate to reality is a gradient distribution for viral particles initially, namely one end of the nasal cavity has more viral particles initially.
In order to take this into account, we simulate the cases where the initial distribution of viral particles follows a gradient. Specifically, the total initial viral particles is still set to be 100. For each row 10 viral particles are assigned. 2 viral particles are assigned to the leftmost 2 grids respectively, 1 viral particle is assigned to the middle 4 grids, and no viral particle is assigned to remaining grids. Similarly, we look at the difference in peak viral load between cases with and without PUFA introduced.
As we can see from figures 24, 25 and 26, the introduction of PUFA results roughly two orders of magnitude lower peak viral load for each grid. Overall, it is safe to conclude that PUFA has made a significant antiviral effect.
.
Conclusions and Future Directions
Using the current parameters sourced from the literature, we initially found that our probiotic is not feasible. This is due to the high clearance rates of our bacterial probiotic and polyunsaturated fatty acids that are in the nasal cavity. We found that adjustment of these parameters yielded a significant effect in terms of our probiotic's ability to decrease viral load. This effect was also found to be robust even after increasing initial viral load, incorporating stochasticity, and taking spatial factors into account. Thus, by tackling the high turnover rate in the nasal cavity, we demonstrated that the smart probiotic could be implemented as an antiviral drug. Our next step was trying to find biologically feasible ways that aim to vary parameters to the desired level.
In order for this probiotic to be feasible, the sensitivity of the synthetic cytokine receptor designed by Aurand and March (2015) would need to be improved. To do this, researchers could use a combination of site-directed mutagenesis and directed evolution. Several mutants of the OprF-OmpA synthetic cytokine receptor can be designed and those with the highest sensitivity and specificity to TNF-α and IFN-γ can be combined, creating a much more sensitive and specific sensor. Protein engineering has been combined with directed evolution to change the specificity and selectivity of other sensors. For example, for the enzyme-based biosensor pyrroloquinoline quinone glucose dehydrogenase (PQQ-GDH), when two mutant forms of the protein were created and combined, the specificity of the sensor became higher for glucose than for other sugars (Yamazaki et al., 2000, Campàs et al., 2009). The same process was done by researchers Dmytruk et al., 2007, in which a mutant form of alcohol oxidase was combined with horseradish peroxidase to create a sensor with a wider linear range than the wildtype form of alcohol oxidase (Dmytruk et al., 2007, Campàs et al., 2009).
In our model, the variations in clearance rate for PUFA was attributed to the incorporation of PUFA into the nasal epithelial cells. Additionally, the clearance rate of PUFA could be overcome by utilizing bioengineering techniques to increase the rates of PUFA production. On our Design page, we discuss techniques that researchers have successfully used to increase de novo PUFA production from synthase complexes.
In addition to increasing PUFA production, the clearance rate of the probiotic cells would need to be decreased. To accomplish this, we consulted ENT Dr. Shikani, who recommended using specific gels in which to deliver the probiotic. For example, researchers could use different types of thermally reversible, mucoadhesive poloxamer gels that have been found to help nasal therapeutics last within the body for longer time periods by adhering slightly to the nasal cavity (Chonkar et al., 2015). These gels contain polymers that increase the retention time in the nasal cavity; they are liquified at cold temperatures and then solidify when placed in the nasal cavity due to warmer temperature acting jointly with the pH level and ions in nasal mucus (Chonkar et al., 2015).
References
Alnaseri, H., Arsic, B., Schneider, J. E., Kaiser, J. C., Scinocca, Z. C., Heinrichs, D. E., & McGavin, M. J. (2015). Inducible Expression of a Resistance-Nodulation-Division-Type Efflux Pump in Staphylococcus aureus Provides Resistance to Linoleic and Arachidonic Acids. Journal of bacteriology, 197(11), 1893–1905. https://doi.org/10.1128/JB.02607-14
Alnaseri, H., Kuiack, R. C., Ferguson, K. A., Schneider, J. E. T., Heinrichs, D. E., & McGavin, M. J. (2019). DNA Binding and Sensor Specificity of FarR, a Novel TetR Family Regulator Required for Induction of the Fatty Acid Efflux Pump FarE in Staphylococcus aureus. Journal of Bacteriology, 201(3), e00602-00618. https://doi.org/10.1128/jb.00602-18
Aurand and March, “Development of a Synthetic Receptor Protein for Sensing Inflammatory Mediators Interferon-g and Tumor Necrosis Factor-a,” Biotechnology and Bioengineering, 2015.
Black P. N. (1988). The fadL gene product of Escherichia coli is an outer membrane protein required for uptake of long-chain fatty acids and involved in sensitivity to bacteriophage T2. Journal of bacteriology, 170(6), 2850–2854. https://doi.org/10.1128/jb.170.6.2850-2854.1988
Campàs, M., Prieto-Simón, B., & Marty, J. L. (2009). A review of the use of genetically engineered enzymes in electrochemical biosensors. Seminars in cell & developmental biology, 20(1), 3–9.
Chonkar, A., Nayak, U., & Udupa, N. (2015). Smart Polymers in Nasal Drug Delivery. Indian journal of pharmaceutical sciences, 77(4), 367–375. https://doi.org/10.4103/0250-474x.164770
Dmytruk, K. V., Smutok, O. V., Ryabova, O. B., Gayda, G. Z., Sibirny, V. A., Schuhmann, W., Gonchar, M. V., & Sibirny, A. A. (2007). Isolation and characterization of mutated alcohol oxidases from the yeast Hansenula polymorpha with decreased affinity toward substrates and their use as selective elements of an amperometric biosensor. BMC biotechnology, 7, 33.
Ericson, M. E., Subramanian, C., Frank, M. W., & Rock, C. O. (2017). Role of Fatty Acid Kinase in Cellular Lipid Homeostasis and SaeRS-Dependent Virulence Factor Expression in Staphylococcus aureus. mBio, 8(4), e00988-00917. https://doi.org/10.1128/mBio.00988-17
Ford, T. J., & Way, J. C. (2015). Enhancement of E. coli acyl-CoA synthetase FadD activity on medium chain fatty acids. PeerJ, 3, e1040. https://doi.org/10.7717/peerj.1040
Green, P. R., Merrill, A. H., Jr, & Bell, R. M. (1981). Membrane phospholipid synthesis in Escherichia coli. Purification, reconstitution, and characterization of sn-glycerol-3-phosphate acyltransferase. The Journal of biological chemistry, 256(21), 11151–11159.
Gullett, J. M., Cuypers, M. G., Frank, M. W., White, S. W., & Rock, C. O. (2019). A fatty acid-binding protein of Streptococcus pneumoniae facilitates the acquisition of host polyunsaturated fatty acids. The Journal of biological chemistry, 294(44), 16416–16428. https://doi.org/10.1074/jbc.RA119.010659
Ishibashi, Y., Aoki, K., Okino, N., Hayashi, M., & Ito, M. (2019). A thraustochytrid-specific lipase/phospholipase with unique positional specificity contributes to microbial competition and fatty acid acquisition from the environment. Scientific reports, 9(1), 16357. https://doi.org/10.1038/s41598-019-52854-7
Kamp, F., & Hamilton, J. A. (2006). How fatty acids of different chain length enter and leave cells by free diffusion. Prostaglandins, leukotrienes, and essential fatty acids, 75(3), 149–159. https://doi.org/10.1016/j.plefa.2006.05.003
Koizumi Y, Iwami S. Mathematical modeling of multi-drugs therapy: a challenge for determining the optimal combinations of antiviral drugs. Theor Biol Med Model. 2014 Sep 25;11:41. doi: 10.1186/1742-4682-11-41. PMID: 25252828; PMCID: PMC4247767.
Lewin, T. M., Wang, P., & Coleman, R. A. (1999). Analysis of amino acid motifs diagnostic for the sn-glycerol-3-phosphate acyltransferase reaction. Biochemistry, 38(18), 5764–5771. https://doi.org/10.1021/bi982805d
Li, G. W., Burkhardt, D., Gross, C., & Weissman, J. S. (2014). Quantifying absolute protein synthesis rates reveals principles underlying allocation of cellular resources. Cell, 157(3), 624–635. https://doi.org/10.1016/j.cell.2014.02.033
Lu, Y. J., Zhang, Y. M., Grimes, K. D., Qi, J., Lee, R. E., & Rock, C. O. (2006). Acyl-phosphates initiate membrane phospholipid synthesis in Gram-positive pathogens. Molecular cell, 23(5), 765–772. https://doi.org/10.1016/j.molcel.2006.06.030
Lu, Y. J., Zhang, F., Grimes, K. D., Lee, R. E., & Rock, C. O. (2007). Topology and active site of PlsY: the bacterial acylphosphate:glycerol-3-phosphate acyltransferase. The Journal of biological chemistry, 282(15), 11339–11346. https://doi.org/10.1074/jbc.M700374200
Lu, X., Vora, H., & Khosla, C. (2008). Overproduction of free fatty acids in E. coli: implications for biodiesel production. Metabolic engineering, 10(6), 333–339. https://doi.org/10.1016/j.ymben.2008.08.006
Maass, S., Sievers, S., Zühlke, D., Kuzinski, J., Sappa, P. K., Muntel, J., Hessling, B., Bernhardt, J., Sietmann, R., Völker, U., Hecker, M., & Becher, D. (2011). Efficient, global-scale quantification of absolute protein amounts by integration of targeted mass spectrometry and two-dimensional gel-based proteomics. Analytical chemistry, 83(7), 2677–2684. https://doi.org/10.1021/ac1031836 Maloy, S. R., Ginsburgh, C. L., Simons, R. W., & Nunn, W. D. (1981). Transport of long and medium chain fatty acids by Escherichia coli K12. The Journal of biological chemistry, 256(8), 3735–3742.
Manganelli, R., & Gennaro, M. L. (2017). Protecting from Envelope Stress: Variations on the Phage-Shock-Protein Theme. Trends in microbiology, 25(3), 205–216. https://doi.org/10.1016/j.tim.2016.10.001
Metz, J. G., Kuner, J., Rosenzweig, B., Lippmeier, J. C., Roessler, P., & Zirkle, R. (2009). Biochemical characterization of polyunsaturated fatty acid synthesis in Schizochytrium: release of the products as free fatty acids. Plant physiology and biochemistry : PPB, 47(6), 472–478. https://doi.org/10.1016/j.plaphy.2009.02.002
Nagano K, Nikaido H. Kinetic behavior of the major multidrug efflux pump AcrB of Escherichia coli. Proc Natl Acad Sci U S A. 2009;106(14):5854-5858. doi:10.1073/pnas.0901695106
Neidhardt and Umbarger, chapter 3 in Neidhardt F.C. Escherichia coli and Salmonella: Cellular and Molecular Biology. 2nd edition. Vol 1. American Society of Microbiology (ASM) Press 1996. table 1
Roberts M. F. (1996). Phospholipases: structural and functional motifs for working at an interface. FASEB journal : official publication of the Federation of American Societies for Experimental Biology, 10(10), 1159–1172. https://doi.org/10.1096/fasebj.10.10.8751718
Rosa, A. O., & Rapoport, S. I. (2009). Intracellular- and extracellular-derived Ca(2+) influence phospholipase A(2)-mediated fatty acid release from brain phospholipids. Biochimica et biophysica acta, 1791(8), 697–705. https://doi.org/10.1016/j.bbalip.2009.03.009
Schönheit, P., Buckel, W., & Martin, W. F. (2016). On the Origin of Heterotrophy. Trends in microbiology, 24(1), 12–25. https://doi.org/10.1016/j.tim.2015.10.003
Schwenk, R. W., Holloway, G. P., Luiken, J. J., Bonen, A., & Glatz, J. F. (2010). Fatty acid transport across the cell membrane: regulation by fatty acid transporters. Prostaglandins, leukotrienes, and essential fatty acids, 82(4-6), 149–154. https://doi.org/10.1016/j.plefa.2010.02.029
Sohlenkamp, C., & Geiger, O. (2016). Bacterial membrane lipids: diversity in structures and pathways. FEMS microbiology reviews, 40(1), 133–159. https://doi.org/10.1093/femsre/fuv008
Stewart, A., Ghosh, M., Spencer, D., Leslie, C. (2002). Enzymatic Properties of Human Cytosolic Phospholiase A2gamma*. The Journal of Biological Chemistry, 227(33), 29526-29536.
Tong, X., Oh, E. K., Lee, B. H., & Lee, J. K. (2019). Production of long-chain free fatty acids from metabolically engineered Rhodobacter sphaeroides heterologously producing periplasmic phospholipase A2 in dodecane-overlaid two-phase culture. Microbial cell factories, 18(1), 20. https://doi.org/10.1186/s12934-019-1070-8
Tonon, T., Qing, R., Harvey, D., Li, Y., Larson, T. R., & Graham, I. A. (2005). Identification of a long-chain polyunsaturated fatty acid acyl-coenzyme A synthetase from the diatom Thalassiosira pseudonana. Plant physiology, 138(1), 402–408. https://doi.org/10.1104/pp.104.054528
Turner, S. M., Moir, J. W., Griffiths, L., Overton, T. W., Smith, H., & Cole, J. A. (2005). Mutational and biochemical analysis of cytochrome c', a nitric oxide-binding lipoprotein important for adaptation of Neisseria gonorrhoeae to oxygen-limited growth. The Biochemical journal, 388(Pt 2), 545–553. https://doi.org/10.1042/BJ20041766
Ujihara, T., Nagano, M., Wada, H., & Mitsuhashi, S. (2014). Identification of a novel type of polyunsaturated fatty acid synthase involved in arachidonic acid biosynthesis. FEBS letters, 588(21), 4032–4036. https://doi.org/10.1016/j.febslet.2014.09.023
van den Berg, B., Black, P. N., Clemons, W. M., Jr, & Rapoport, T. A. (2004). Crystal structure of the long-chain fatty acid transporter FadL. Science (New York, N.Y.), 304(5676), 1506–1509. https://doi.org/10.1126/science.1097524
Wolf, M. J., & Gross, R. W. (1996). Expression, purification, and kinetic characterization of a recombinant 80-kDa intracellular calcium-independent phospholipase A2. The Journal of biological chemistry, 271(48), 30879–30885. https://doi.org/10.1074/jbc.271.48.30879
Wu, T., Ikezono, T., Angus, C. W., & Shelhamer, J. H. (1994). Characterization of the promoter for the human 85 kDa cytosolic phospholipase A2 gene. Nucleic acids research, 22(23), 5093–5098. https://doi.org/10.1093/nar/22.23.5093
Xiao, K., Yue, X. H., Chen, W. C., Zhou, X. R., Wang, L., Xu, L., Huang, F. H., & Wan, X. (2018). Metabolic Engineering for Enhanced Medium Chain Omega Hydroxy Fatty Acid Production in Escherichia coli. Frontiers in microbiology, 9, 139. https://doi.org/10.3389/fmicb.2018.00139
Yamazaki, T., Kojima, K., & Sode, K. (2000). Extended-range glucose sensor employing engineered glucose dehydrogenases. Analytical chemistry, 72(19), 4689–4693.
Zitzmann, C., & Kaderali, L. (2018). Mathematical Analysis of Viral Replication Dynamics and Antiviral Treatment Strategies: From Basic Models to Age-Based Multi-Scale Modeling. Frontiers in microbiology, 9, 1546. https://doi.org/10.3389/fmicb.2018.01546