Frequency and Time Domain Analysis ofsRNA-based Treatment for Inflammatory BowelDisease
Arman Karshenas, Joseph Windo, Bhuvana, Jonathan Stocks, Adrian, Kozhevnikov, Jhanna Kryukova, Eleanor Beard, Laurel Constanti Crosby, Max, Taylor

TL;DR
This paper validates a complex reaction pathway for treating Inflammatory Bowel Disease using frequency and time domain analysis, highlighting the potential of sRNA-based probiotics for future therapeutic applications.
Contribution
It introduces a novel frequency and time domain analysis approach to verify a proposed reaction pathway for IBD treatment, emphasizing the role of negative feedback loops.
Findings
Negative feedback loop significantly influences the reaction pathway.
Proposed probiotics show promising potential for IBD treatment.
Frequency domain analysis confirms pathway validity.
Abstract
The validity of a complex reaction pathway proposed to treat Inflammatory Bowel Disease (IBD) was verified by a comprehensive time and frequency domain analysis. The model was taken to the frequency domain to study the effect and the significance of the negative feedback loop introduced by the reaction pathways. It could be shown that such proposed probiotics have very interesting potentials that could be used extensively in the near future.
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
TopicsAdvanced biosensing and bioanalysis techniques · RNA and protein synthesis mechanisms
Frequency and Time Domain Analysis of sRNA-based Treatment for Inflammatory Bowel Disease*
††thanks: BBSRC, Wellcome Trust, IDT, Eppendorf, BioLab, Snap gene and Oxford university departments of biochemistry, Engineering and chemistry
1st Arman Karshenas
Department of Engineering Science
*University of Oxford
*Oxford, UK
2nd Joseph Windo
Department of Biochemistry
*University of Oxford
*Oxford, UK
3rd Bhuvana Sudarshan
Department of Biochemistry
*University of Oxford
*Oxford, UK
4th Jonathan Stocks
Deptartment of Biochemistry
*University of Oxford
*Oxford, UK
5th adrian Kozhevnikov
Department of Engineering Science
*University of Oxford
*Oxford, UK
6th Jhanna Kryukova
Department of Biochemistry
*University of Oxford
*Oxford, UK
7th Eleanor Beard
Department of Medicine
*University of Oxford
*Oxford, UK
8th Laurel Constanti Crosby
Department of Biology
*University of Oxford
*Oxford, UK
9th Max Taylor
Department of Biochemistry
*University of Oxford
*Oxford, UK
Abstract
The validity of a complex reaction pathway proposed to treat Inflammatory Bowel Disease (IBD) was verified by a comprehensive time and frequency domain analysis. The model was taken to the frequency domain to study the effect and the significance of the negative feedback loop introduced by the reaction pathways. It could be shown that such proposed probiotics have very interesting potentials that could be used extensively in near future.
Index Terms:
Synthetic Biology, Control theory, Negative feedback in Biological systems
I Introduction
Inflammatory Bowel Disease (IBD) is characterised by chronic inflammation of the intestine. The condition is associated with an imbalance in immune cell populations, notably Th17 and Treg. Existing immunosuppressive therapies, when successful, often elicit systemic side effects and require frequent re-administration. Therefore, the proposed solution is a probiotic strain that restores the Th17/Treg cell balance via secretion of IL-10 in response to Nitric oxide in the intestinal lumen. Overshoot is prevented by an adenine riboswitch-sRNA construct which responds to extracellular adenosine, an indicator of the Treg cell population. The two independent reaction pathways - responding separately to IL-10 deficiency and abundance - and are summarised in Figure 1.
To avoid the complexity of multivariate control systems, the pathways have been modelled independently. Physiologically, an elevated serum concentration of NO can be associated with autoimmunity, whereas an elevated serum adenosine concentration can signify immunoinsufficiency. The NO pathway will be referred to as the ‘negative pathway’, representing a negative concentration difference in IL-10. Likewise, the adenosine pathway has been termed the ‘positive pathway’. Integration of separate stimuli in a dual feedback loop enables a more dynamic, robust response to the immune state of the body. Thus the engineered probiotic bacteria would be capable of maintaining an equilibrium between the Th17 and Treg cell populations by virtue of the negative feedback system.
I-A Negative Pathway
The negative pathway can be modelled in a canonical manner: nitric oxide triggers oxidation of the SoxR iron-sulphur cluster, allowing binding to the promoter of IL-10 in the circuit - pSoxS. This results in transcriptional activation of the system within approximately seconds [1]. Hence, the system has been modelled by a simple Hill activation function for IL-10 gene transcription. The SoxR/pSoxS system has been studied extensively [2] but previous modelling has been limited to that of the SoxR/pSoxS binding interaction by the 2009 Stanford iGEM team [3].
I-B Positive Pathway
The positive pathway mechanism facilitates translational inhibition of IL-10 to prevent overshoot. The anti-inflammatory molecule, adenosine, is produced by Treg cells in response to inflammation, with elevated levels indicating Treg overactivity and an immunoinsufficient state. Adenosine is subsequently broken down into ribose and adenine, via an extracytoplasmic, nucleoside hydrolase, by the engineered cell. Following uptake into the cell, adenine binds to a riboswitch, resulting in an increased rate of RNA transcription. A self-splicing ribozyme separates the riboswitch from the remaining transcript - an sRNA molecule which is complementary to the 5’ region of the IL-10 mRNA.
Consequently, translation of IL-10 mRNA is prevented by the binding of sRNA. It should be noted that NO is an essential input to the positive pathway and hence, combined model for both reaction pathways could be found below in Figure‘2.
II Dynamic Model
The dynamic model of both of the pathways were analysed using the Simbiology toolbox in MATLAB as well as ODE15s Solver and some of the modelling and parameters are based on the paper ”Frequency domain analysis of small non-coding RNAs” [4]. The dynamic time domain analysis of the system was conducted using the combined model in Figure 2 with different initial conditions corresponding to negative and positive pathways separately. It is worth mentioning that data fitting was used to determine some of the parameters. A list of the parameters and data is available at the end of paper. The differential equations describing the systems are as follows:
[TABLE]
[TABLE]
[TABLE]
[TABLE]
where are the Hill functions describing the rate of transcription with inputs of Adenosine and NO for and respectively as represented by in (5).
[TABLE]
In which is the maximal transcription factor, the dissociation coefficient and the Hill coefficient. is the rate of sRNA binding to IL-10 mRNA which heavily depends on the length of the sRNA, making it easy to exploit this relationship for better control over the fate of the system. is the translation rate of IL-10 mRNA and is the secretion and diffusion rate of IL-10 [5]. represents the dilution + degradation rate.
It is important to notice that Hill functions introduce non-linearity into the system as well as sRNA binding. Therefore, Hill functions have been linearised for transfer function derivation in the frequency domain, and are noted by as follows:
[TABLE]
All the assumptions made for (1) to (4) are listed below:
- •
Adenosine substituted with adenine as the hydrolase reaction is believed to be much faster than the body response
- •
Concentration of adenine and NO kept constant for dynamic analysis due to intra- and extracellular abundance
- •
The initial conditions used for time domain analysis correspond to maximal insufficiency based on [6], [7] and [8].
- •
Stochastic response ignored due to the large number of E. Coli used.
A list of all parameters and concentrations used in figures, charts and graphs can be found at the end of the section.
II-A Negative pathway dynamics
Negative pathway dynamics are fully described by the ODEs stated previously. It should be noted that based on the second assumption, the level of NO would not change with time and hence . An initial concentration of was used for elevated level of NO corresponding to IL-10 deficiency in patients with IBD [6] as well as the nominal concentration of for adenine [8]. The evolution of the system to NO stimuli is illustrated in Figure 3.
As it can be seen from Figure 3, reaches a concentration of after roughly 180 seconds. It is important to note that the initial concentration of IL-10 is set to 0 for simplicity, meaning that in reality, the model can correct deficiencies of IL-10 to a margin of . A control mechanism is essential in order for the probiotic bacteria to be able to maintain a healthy level of IL-10. To determine how system would respond and whether or not it would be able to stabilise the IL-10 concentration, it is necessary to construct a large-scale version of the model that incorporates the body as well as the system embedded in a negative feedback loop [9]. The control modelling is discussed in more detail later in the paper. Lastly, the model behaviour is as expected due to the constant production rate, as opposed to the degradation rate, and hence a bounded steady state value is expected.
II-B Positive pathway dynamics
The complexity of the positive pathway is visible from Figure 2. The complexity of the system was maintained for the dynamic analysis based on ODEs stated previously. It should be noted that based on the second assumption, the concentration of adenine would not change with time and hence . An initial concentration of was used for elevated level of adenine corresponding to IL-10 abundance in patients with IBD [8] as well as the nominal concentration of for NO [6]. The evolution of the system to adenine stimuli is illustrated in Figure 4.
It can be observed from Figure 4, that is inhibited significantly and it could be assumed that a high concentration of adenine could inhibit translation of IL-10 completely. It should be noted that the length of sRNA and the promoter strength in both negative and positive pathways() play a significant role and hence they have been optimised for maximum inhibition as well as maximum translation of IL-10 in both positive and negative pathways respectively. The sensitivity analysis and sRNA models are discussed below.
III Sensitivity analysis and sRNA modelling
Sensitivity analysis was done separately for each reaction pathway in order to find the most significant steps and parameters, thus enabling these parameters to be exploited in order to improve the responsiveness of the system to smaller perturbations in concentrations. The sensitivity of with respect to parameters used in simulations is shown in Figure 5.
It is clear from Figure 5 that IL-10 production depends heavily on and which is desirable as different promoter strengths (weak, medium and strong) would enable optimisation of the response. It can also be observed that sRNA binding rate has less significant effect on the output compared to and hence, a stronger promoter would mean that a much longer sRNA would be required to inhibit IL-10 translation. The length of sRNA was optimised to 24 base pairs [10] yielding . The graph of against the position of base pairs is illustrated in Figure 6.
Based on the calculated , , and were set to , and for complete inhibition of IL-10 translation in the positive pathway, in addition to a high IL-10 translation rate in negative pathway.
IV Steady state model
The simplest model is an output/input, time-independent plot which describes the system fate in response to varying input concentrations. Secretion of IL-10 has been ignored for steady state response as it is assumed that the system is given ample time to settle, secreting all translated IL-10. The steady state concentrations of species can be determined by equating (1), (2) and (3) to zero. The steady state values for [IL-10], [sRNA] and [IL-10 mRNA] are given below:
[TABLE]
[TABLE]
[TABLE]
where \omega=\big{(}k_{1}(\Gamma_{2}-\Gamma_{1})-\alpha_{1}\alpha_{2}\big{)}^{2}+4\alpha_{1}\alpha_{2}k_{1}\Gamma_{2}. It should be noted that denotes Hill functions for adenine and NO. A surface can be plotted describing the output/input steady state relation, with contours representing the steady state characteristic of each of the pathways. The plots can then be used for experimental fitting to obtain more realistic parameters. The steady state surface and it’s contours could be found in Figure 7.
V Parameters and modelling predictions for circuit design
The parameters used for modelling are all listed below:
- •
Transcription level : , both optimised to give satisfactory results [4]. [11], [4].
- •
Translation level : optimised based on the number of base pairs in sRNA. is the translation rate of IL-10 [4].
- •
Degradation and Diffusion : , , all based on [4]. The diffusion rate parameter was calculated based on the number of transport proteins on the membrane of E. Coli to be [12].
- •
Body response : The value of and were calculated from Figure 8 and a chosen value of s. and .
V-A Parameters range and output dependency
Three of the parameters used throughout the paper can be tuned for optimisation, namely: , and ; some of which have been previously mentioned. These parameters were modified slightly about their optimised value to find a range for which the system maintains a desired working state. The optimised range for each of the parameters is given below.
has a wide range of for which a satisfactory results could be obtained. 2. 2.
As opposed to , has a very narrow range of for which a satisfactory results could be obtained. 3. 3.
is almost saturated at the value chosen for simulation and hence does not have an upper bound; however, there is a lower bound of for which a satisfactory result can be obtained.
Equations 7 - 9 can be differentiated to find the dependency of IL-10 with respect to all three tunable parameters. As depends linearly on , (7) has been differentiated with respect to , and for simplicity. It should be noted that since is the only tunable parameter in Hill functions, has been used as the parameter for simplicity in differentiation.
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
where \lambda=\frac{\big{(}k_{1}(\Gamma_{1}-\Gamma_{2})-\alpha_{1}\alpha_{2}\big{)}(\Gamma_{1}-\Gamma_{2}\big{)}+4\alpha_{1}\alpha_{2}\Gamma_{2}}{\sqrt{\omega}} and \mu=\big{(}\alpha_{1}\alpha_{2}+k_{1}(\Gamma_{1}-\Gamma_{2})+\sqrt{\omega}\big{)}. Based on (10), (11) and (12), the dependency of IL-10 with respect to , and is clear. The expression is less than and therefore it could be concluded that:
2. 2.
3. 3.
and are both in the inhibitory pathway as opposed to and hence, all three statements above validate the model.
V-B Modelling predictions for circuit design
Based on the range and dependency of on different parameters, a medium strength promoter is suggested to be used for negative pathway as opposed to positive pathway where either a strong promoter or a medium strength promoter with a high copy number is suggested. As for the length of sRNA, it was shown that as long as it has just a few base pairs, corresponding to , the output would not be affected significantly and hence, the most energetically feasible length is used.
VI Body Dynamics
The response of the body is crucial to obtain comprehensive model reflecting systemic interactions between the engineered E. coli and the body. Body response mechanisms have been studied in detail however due to stochastic noise, significant assumptions were made to simplify the response for modelling. These assumptions are summarised below; please note that even though some of the assumptions are considerable, they are negligible compared to the amount of variation present in bodily systems and between individuals:
- •
The model is based on bulk amount of bacteria translating or inhibiting IL-10.
- •
The body responds in about 3-4 hours and the response has been taken to be linear due to lack of data points.
- •
The body has been assumed to behave similarly to an LTI (Linear time-invariant) system and the model can be extrapolated, with care, for people with different body responses.
Based on assumptions stated above, two lines can be plotted for elevated levels of NO and Adenosine versus IL-10 - the highest and lowest concentration of IL-10 [13] and its corresponding signalling molecule [6] [8] that has been measured in human body. The plot of steady-state concentrations is shown in Figure 8.
The dynamic response of the body is assumed to be linear with a response time of about 3 hours as very few data points were available for the steady state response of the body. Therefore, the differential equation describing the response is as follows:
[TABLE]
Figure 8 can be used to find and for both positive and negative pathways. Then by choosing , for the response time, (13) is fully defined and can be solved simultaneously with (1) - (4) to give an overall prediction on how the body and the probiotic bacteria would interact with each other to reach an equilibrium concentration of IL-10 in the body. It should be noted that for these two to be solved together, one should be aware of the correction factor required to approximate metabolite concentrations inside a bacterium relative to the body, in addition to an indication of the number of bacteria being used. Hence, the initial number of bacteria is optimised in a correction loop by looking at the settle time of IL-10 inside the body. The dynamics of IL-10 in the body for elevated level of NO and adenine are shown in Figures 9 and 10 respectively. It should be noted that for representation purposes, the concentration of IL-10 in the body - magenta dashed line - has been multiplied by .
It is clear from Figures above that the body and E. Coli reach a stable concentration of IL-10 without oscillation in any of the metabolite concentrations. The accuracy of the response time of seconds is reflected by the similarity of the steady-state.
The number of E Coli used for simulation would significantly affect the settle time or the response time and hence, the time taken by the combined system to reach of the steady-state concentration of IL-10 has been plotted against the number of bacteria in Figure 11. It is essential to understand the significance and importance of the positive pathway in stability and robustness of the proposed probiotic. Hence, dynamics of the combined system has been plotted in absence of Adenine pathway to demonstrate that the system becomes unstable without a negative feedback loop.
VII Negative feedback and frequency domain analysis
As discussed earlier in the paper, both positive and negative feedback could potentially correct the concentration of IL-10 in the body. However, it is important to note that body response should be included in the model to give an actual and meaningful prediction of how E. coli would behave in the body. For this purpose, control theory and the negative feedback loop have been used by transforming the body and the bacteria in a standard cascade compensation model and deriving their transfer functions. Hence, we developed a cascade compensation model where we have the body as the plant and the bacteria as the controller. Assumptions made for frequency domain analysis are listed below:
- •
The system has been linearised about an equilibrium point and therefore, it has been assumed that perturbations are small for the Taylor series to converge.
- •
The reaction pathways will be analysed separately as based on superposition, the response of the linearised model would be a sum of its response to inputs separately.
- •
Secretion has been ignored for simplicity in transfer function.
Therefore, by using the correct transfer function, both models would fall into the same control design, just with different transfer functions. The proposed design for the feedback loop could be seen in Figure 13.
The aim of the controller is to set the reference signal, namely NO and adenosine, to nominal level as this would correspond to a healthy concentration of IL-10 within the body. Transfer functions for the controller- both positive and negative pathway- and the body are required for control analysis of system response in a human body.
VIII Transfer functions
As stated in the assumptions, the model has been linearised due to non-linearity raised by the Hill functions and sRNA binding. Therefore, Hill functions have been replaced by as stated earlier in Equation 6. The block diagram and signal flow graph for the linearised system are shown in Figures 14 and 15, where transcription of IL-10 mRNA, sRNA and translation of IL-10 can be seen clearly.
Transfer functions for both negative and positive pathways were found using Mason’s gain relation [14] and confirmed against [4]. Transfer functions of both positive and negative pathways can be found below:
[TABLE]
[TABLE]
where and .
VIII-A Body response transfer function
Laplace transform can be taken from both sides of the Equation 13 and by ignoring , being very small and negligible, the transfer function for the body can be derived:
[TABLE]
[TABLE]
IX Stability
The closed loop transfer function for cascade cascade controller design given in Figure 13 would contain important information regarding steady state error and relative stability. It should be mentioned that similar to most control problems where the controller is designed for specific specifications [15], the reaction pathway can be changed by altering promoter strength and sRNA length, making it almost possible to adopt desired structures for the controller. The closed loop transfer function for positive and negative pathways is as follows:
[TABLE]
Based on (17), the stability of the system can be determined by observing the location of the poles which are roots of the characteristic polynomial . The characteristic equation for adenine and NO reference signals is given in Equations 18 and 19.
[TABLE]
[TABLE]
Equations (18) and (19) were solved numerically [16] and smallest poles and zeros for both reaction pathways are summarised in the Argand diagram below, Figure 16. All poles and Zeros are summarised after the figure.
It should be noted that some poles and zeros were not plotted as they were significantly larger(in magnitude) than those shown in Figure 16 and would have thus hindered comparison with those illustrated. All poles and zeros are given below.
NO response :
- •
Poles : , , and .
- •
Zeros : , , , . 2. 2.
Adenine response :
- •
Poles : , , and .
- •
Zeros : , , , .
Despite right-plane poles, it should be noted that both responses can be considered stable as the dominant poles are significantly larger- in magnitude- than all other poles.
X Steady state error - DC gain
[TABLE]
[TABLE]
[TABLE]
By substituting back Transfer functions for both E. Coli,, and the body,, the steady state error can be found as follows:
- •
For NO response .
- •
For adenine response .
It is evident that the combined model has a very small DC gain and hence, a very small margin of error is estimated for the steady concentration of IL-10.
Acknowledgement
Supervision and guidance of the following people were significant in achieving and finishing this paper. Their contributions and help is greatly appreciated.
- •
Karandip Saini - Department of Chemistry, our team wet Lab coordinator
- •
Dr. Aivar Sootla - Department of Engineering sciences University of Oxford.
- •
Dr.Nicolas Delalez - Department of Engineering sciences University of Oxford.
- •
Prof.Antonis Papachristodoulou - Department of Engineering sciences University of Oxford.
- •
Prof.Christopher Macminn - Department of Engineering sciences University of Oxford.
- •
Dr.George Wadhams - Department of Biochemistry University of Oxford.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1[1] Uri Alon. An Introduction to System Biology . Chapman, 2006.
- 2[2] Bruce Demple. Redox signaling and gene control in the escherichia coli soxrs oxidative stress regulon- a review. GENE , 179:53–57, 1996.
- 3[3] Stanford igem 2009.
- 4[4] Harrison Steel and Antonis Papachristodoulou. Frequency domain analysis of small non-coding rnas shows summing junction-like behaviour. IEEE , (56):1014–1023, 2017.
- 5[5] Michael H. H. Lenders, Tobias Beer, Sander H. J. Smits, and Lutz Schmitt. In vivo quantification of the secretion rates of the hemolysin a type i secretion system. Nature scientific reports , (6), 2016.
- 6[6] Nesina Avdagic and Asija Zaciragic. Nitric oxide as a potential biomarker in inflammatory bowel disease. Bosnian Journal of basic medical sciences , (13):5–9, 2013.
- 7[7] Keiichi Mitsuyama, Nobuo Tomiyasu, and Kosuke Takaki. Interleukin-10 in the pathophysiology of inflammatory bowel disease: Increased serum concentrations during the recovery phase. Mediators of inflammation , 2006.
- 8[8] Thomas W. Traut. Physiological concentrations of purines and pyrimidines. Molecular and Cellular Biochemistry , (140):1–22, 2006.
