Model
"The essence of mathematics lies in its freedom."
View ENTIRE REPORT of our model which embraces all analysis details and derivation processes of the formulas and calculation methods.
Use our SOFTWARE VisualSTEP to explore the model freely!
The core of our project is the signal transducer system on the cell membrane, which consists of a transcription factor chain (TC chain) and a protease (PR) chain (PC chain). The two chains can be combined with the help of an exogenous ligand, thus producing an transcription factor (TF), which is regarded as the output. This is the main mechanism of our system[1], as is shown in Figure 1.
There are several properties measuring the faculty of the system: background noise, system (saturated) output, dynamic range (DR) and system efficiency.
The system output without exogenous ligand, regarded as the background noise, comes from the random collision between the two chains[1]. And the system saturated output, on the contrary, characterizes the maximum response strength of our system. Without question, the system will be superior if we maintain the saturated output as well as reduce background noise, namely the increasing of dynamic range. Moreover, the TC chain without transcription factor may still be combined with a PC chain, which makes the PC chain occupied but does not contribute to system output. This phenomenon reduces the utilization of system resources, which is measured by the indicator called system efficiency. The practical significance of promoting system efficiency lies in achieving the same level of system output with fewer resources (i.e., lower chain generation rates).
Figure 1. The mechanism of our system.
Based on the above, this is the main content of our modeling:
• We establish a differential equation model based on the system mechanism described above. As is partly verified by our wet lab, we quantitatively prove that the maximum of dynamic range can be reached by changing the production rates of PC and TC chains[1] and demonstrate how to do it.
• For the high background noise from wet lab results, the model above successfully explains why the excessively high chain production rates enhance background noise. Furthermore, we provide a method, that can be used in further experiment, to simultaneously achieve greater saturated output and lower background noise.
• We also establish a model based on birth-death process, so that we quantitatively prove that the dynamic range is negatively correlated with protease cleavage ability under normal condition, which gives a explicit guidance on further experiment.
• Furthermore, We advance the birth-death process to a more general continuous-time Markov chain, which demonstrate that it is feasible to increase system efficiency by changing binding abilities, as long as we design a dedicated degradation pathway of free TC chains without transcription factor. What’s more,we obtain a method to improve system efficiency while preventing system saturation output from decreasing. This will be one of the ideas for our further experiment.
In the wet experiment, we can control the expression rates of TC chain and PC chain by regulating the content of plasmids transferred into the cells, and then indirectly control the quantitative ratio of the two chains on the membrane. Although the expression rate and the number of receptor chains are not strictly proportional, there is no doubt that there is a positive correlation between them.
We establish a system of differential equations that can reflect the reaction mechanism of the entire system in detail and use it to quantitatively prove that the dynamic range can be improved by changing the production rates of PC and TC chains, under normal condition. In addition, via the model, we provide a specific method of changing the two values mentioned above.
In the first place, we assume the operation of the system can be described by a seire of chemical equations, and obtain the corresponding differential equations, through law of mass action. The production rates of the two chains can be balanced with a dilution rate that all substances in the system have, so the system can reach a stable state after a while. The initial differential equations and a simulation with appropriate parameters are both shown in Figure 2.
There are product terms of independent variables in the initial equation set so that we can only obtain the numerical solution of it. We use the linear approximation system, obtained by the Taylor series expansion of those product terms, to find the analytical solution which is necessary for our exploring the relationships between system properties and parameters. Although Lyapunov Stability Theory states that the stability of the primary system can be represented by the stability of its linear approximation system only if some conditions are met[6][7], there is no doubt that the whole system will eventually reach a steady state[4].
By letting the equations of the linear approximation system to be equal to zero, we get the function about the relationship between dynamic range and the production abilities of the two chains. In order to be more realistic, we consider the effects of random errors, downstream leakage measurement errors and other factors on dynamic range, which are not significantly related to experimental parameters. In ohter words, we can regard their mean values as constants and consequently get Figure 3, in which we show the original dynamic range function on the left and the graph of ultimate DR function on the right.
Figure 2. The initial differential equation set and its simulation. For convenience, we use A to represent for the PC chain, B to represent for the TC chain carrying the transcription factor, P to represent for their coalition, and R to represent for the transcription factor. Meanwhile, we use the footnote "l" to indicate that the substance has lost the transcription factor. That is, a TC chain out of transcription factor represented by Bl, and a coalition of Bl and A is represented by Pl.
Figure 3. The DR function and its graph. The horizontal and vertical coordinates indicate the production rate of PC and TC chains, respectively(the transfection concentrations of the two chains can be used to replace them in wet experiment), and the color indicates the size of DR.
The figure shows the existence of the maximun value of dynamic range, which is verified by the wet experiment. For instance, we can get the curve of DR when the ratio of the horizontal and vertical coordinates is constant, which is equivalent to the vertical cutting through the origin of Figure 3, as shown in Figure 4.
In conclusion, our model reveals the intrinsic link between production rates and dynamic range. Thanks to it, we now have a clearer way to achieve a higher dynamic range.
Figure 4. A way to maximize dynamic range.
In the beginning, our wet experiment shows a relatively large output of our system in spite of low concentration of the exogenous ligand, which is not quite as expected. This phenomenon reduces the dynamic range of the system, resulting in a decrease in the sensitivity of the system to the concentration of exogenous ligands.
Through the formula we obtained from the modeling of differential equations (Equation 18 and 19 in Section 1.2.2 of our entire model report), α and β are both positively correlated with the system output. This explains why the excessively high chain generation rate enhances the background noise. So we can simultaneously achieve greater dynamic range and smaller background noise by reducing the two production rates at the same time, as long as the factor (α-a0ζ)/β is always small.
Due to random collision or exogenous ligand, a PC chain and a TC chain will be combined. There are two situations after that. In the first case, after a period of time T the protease cleaves the transcription factor, then the two chains are separated. In the second case, there are not enough time for the protease to cleave the transcription factor before the two chains are separated. Obviously T is a random variable, and the more likely the cleavage occurs, the smaller its mathematical expectation E(T) is.
Whether the cleavage is prone to occur is mainly determined by the catalytic efficiency of the protease, and is also affected by other factors such as the length of the intracellular segment of the PC chain and the TC chain. For the sake of convenience, we say that protease cleavage ability is stronger if the cutting is more likely to occur.
In order to research the relationship between dynamic range and protease cleavage ability via this stochastic process, we need to consider two typical cases: the concentration of exogenous ligand is nearly 0 or large enough so that the system output is saturated.
Let us focus on how a TC chain behaves in the two situations. At first, it moves freely on the membrane (State 1, the free sate). After a period of time it is combined with the PC chain (State 2, the bound state) due to random collision (in the first case) or an exogenous ligand (in the second case), and then it is either cleaved by protease (State 3, the inactive state) or separated from the PC chain without being cut (returns to State 1). Once it reaches State 3, it loses the transcription factor, so whether it is combined with the PC chain, it will no longer increase the output of the system. In other words, it cannot transition from state 3 to state 1 or 2.
Based on the above content, we establish a birth-death process model, whose state transition diagram is shown in Figure 5.
Figure 5. State transition diagram of the birth-death process. The TC chain can transition its states along the arrows with the corresponding exponential distribution parameters p,c and d. Since the mean of a exponentially distributed variable is the reciprocal of its exponential distribution parameter, the longer the mean value is, the smaller the corresponding probability parameter is. Therefore, in the case of stable external environment (temperature, pH, etc.), p is positively correlated with the concentration of PC chain, d is negatively correlated with the binding ability of the two chains, and c is positively correlated with protease cleavage ability.
The mean of the total time E(Me) that a TC chain has undergone from being placed on the cell membrane to being cleaved by the protease reflects the activity of the system (i.e., the greater the output of the system is, the bigger the reciprocal of it is). Hence, dynamic range is positively correlated with the expression shown in Figure 6. We also display the graph of it with protease cleavage ability as the independent variable.
It is clear from the figure that the negative correlation between dynamic range and protease cleavage ability is caused by different growth rates of background noise and saturated output when protease cleavage ability rises, which gives a explicit guidance on further experiment.
Figure 6. The relationship between dynamic range and protease cleavage ability. In order to show the relationship between protease cleavage ability and dynamic range better, we assume DR=kE(Me1)/E(Me2) without loss of generality, where k is a scale factor.
After losing its transcription factor, a TC chain may still be combined with a PC chain, which makes the latter occupied but does not contribute to system output.
We use the indicator called system efficiency (SE) to measure the extent, to which TC chains out of transcription factor occupy system resources.
By promoting SE, we can achieve the same level of system output with fewer resources (i.e., lower chain generation rates).
It should be noted that since the TC chain is no longer likely to contribute to system output after losing its transcription factor, using the TC chain to measure SE is a bit of a hassle.
However, the structure of the PC chain is always the same, which means that we can easily measure $SE$ with parameters related to PC chain.
Therefore, we advance the birth-death process model above to an ordinary continuous-time Markov chain, in order to research the behavior of a PC chain. The state transition diagram of the advanced model is shown in Figure 7. The meaning of each state the PC chain can be in is as follows:
State 1: (free state) it moves freely on the cell membrane,
State 2: (semi-bonded state) it binds to an exogenous ligand,
State 3: (effective binding state) it binds to a TC chain with the help of exogenous ligand, so that they may produce transcription factors,
State 4: (invalid binding state) it binds to a TC chain without transcription factor with the help of exogenous ligand, which means it is occupied but cannot help produce transcription factors.
Using the theory of stochastic process[5], we find the analytical solution of the stability probability of this chain in each state after a long enough time, which can help measure system efficiency. Hence, we obtain the relationship between the binding abilities of the two chains to the exogenous ligand.
Figure 7. State transition diagram of the Markov process.
Through the equations obtained from the model (Eqs. 42, 43 and 44 in Section 2.2 of our entire model report), it follows that under the condition that our basic system does not contain a specific degradation pathway to free TC chain out of transcription factor, we cannot increase system efficiency by changing binding abilities. However, if we design the pathway, it will be feasible to increase system efficiency by improving the binding ability of the PC chain and exogenous ligand or reducing that of the TC chain and exogenous ligand, which is shown in Figure 8.
Moreover, improving the binding ability of the PC chain and exogenous ligand can help promote the combination of TC and PC chains, which leads to a higher output. It follows that we obtain a method to improve the system efficiency while preventing system saturation output from decreasing: improving the binding ability of the PC chain and the exogenous ligand while reducing the binding ability of the TC chain and the exogenous ligand at the same time. This will be one of the ideas for our further experiment.
Figure 8. Effect of the binding abilities on system efficiency under the condition of the existence of a degradation pathway to free TC chain out of transcription factor. λ1 are negatively correlated with the binding ability of the PC chain and exogenous ligand, and λ2 are negatively correlated with that of the TC chain and exogenous ligand.
View ENTIRE REPORT of our model which embraces all analysis details and derivation processes of the formulas and calculation methods.
Use our SOFTWARE VisualSTEP to explore the model freely!
[1]Schwarz, K. A., et al. Rewiring human cellular input-output using modular extracellular sensors. Nature Chemical Biology 13.2(2016):202.
[2]Daringer, N. M., et al. Modular extracellular sensor architecture for engineering mammalian cell-based devices. Acs Synthetic Biology 3.12(2014):892-902.
[3]Schwarz, Kelly A., and J. N. Leonard. Engineering cell-based therapies to interface robustly with host physiology. Advanced Drug Delivery Reviews 105.Pt A(2016):55-65.
[4]Oppenheim, Alan V., et al. Signals and Systems, 2nd Edition. Englewood: Prentice Hall, 1996.
[5]He, Shuyuan, et al. Stochastic Process. Beijing: Peking University Press, 2008.
[6]Cai, Suilin, et al. Ordinary Differential Equation, 3rd Edition. Hangzhou: Zhejiang University Press, 2012.
[7]Liu, Wending, et al. The Principle of Automatic Control, 3rd Edition. Beijing: Publishing House of Electronics Industry, 2013.
Address
G604, School of Life Sciences, Fudan University
2005 Songhu Road, Yangpu, Shanghai, China
2005 Songhu Road, Yangpu, Shanghai, China