Team:McGill/Modeling

From 2009.igem.org

Latest revision as of 01:54, 22 October 2009

The purpose of the following mathematical modeling is to gain insight into activation-inhibition signaling. We begin by describing the mathematical model employed to simulate our biological system. Then we examine the dynamics of one activation site and one inhibitory site, what we call an oscillator, and the effect of distance. We also examine the role two parameters, the Hill coefficient and diffusion rate, have on the oscillatory behavior. This leads to the next examination of two oscillators and the effect of varying the distance between the oscillators but holding the distance between the sites within an oscillator constant. Here we find several unexpected dynamics that we further explore. Finally we look into two dimensional modeling. However because of the size of the .avi files we will not discuss this modeling here but save it for our presentation at the Jamboree.

The following represents a brief overview of all of our modeling exercises. We unfortunately did not make use of this wiki throughout the summer and so many of our results have been left out for lack of time to upload. However, all of our results will be available on our poster and oral presentations at the Jamboree.

Introduction

Many models examining intercellular signaling do not take into account the separation distances of the signaling bodies. We use a partial differential equation (PDE) based model to gain insight into spatially heterogeneous activation-inhibition intercellular signaling.

Two types of signaling molecules exist: activating and inhibiting. Each molecule is synthesized by a unique strain of cells and affects the synthesis rate of the other strain.

Figure 1 – Activation-inhibition intercellular signaling – Activating molecule (A) synthesizes and diffuses to increase synthesis of inhibiting molecule (B) in secondary strain. Inhibiting molecule also diffuses back to initial cell and decreases synthesis of activating molecule.

This is modeled using the following system of PDEs:

where Ψ1 and Ψ2 represent the concentrations of the activating and inhibiting molecules, respectively, γi the degradation constant, Di the diffusion constant, λi the maximal synthesis rate of molecule i, and δ the Dirac function. fi represents the Hill function describing the dependence on the opposing molecule:

where n, b, and θ are positive. To simplify the analysis, we chose equal parameters between the activating and inhibiting sites (Appendix A).

Numerical Simulation

The above system was solved numerically using a forward Euler scheme in time and a centered difference scheme in space. Cyclical boundary conditions were assumed; meaning the spatial dimension formed a ring. This was chosen since simulating an approximate infinite line is computationally costly. However, this also allowed us to investigate two site geometries during one simulation (explained later). The ring was given a physical length of 50 and discretized into 500 intervals. For simplicity, separation distances will be reported in terms of numerical intervals rather than physical distance.

One Oscillator

We first explored the different potential dynamics when the separation distance between an activating and inhibitory site was increased.

Figure 2 – One Oscillator – The red bar represents the inhibitory site, which remains fixed in position while the activating site, blue bar, is sequentially moved around the ring.

The following is an example of the dynamics observed when the two sites are at a distance where oscillations occur.

Figure 3 – One Oscillator – One activation and one inhibitory site separated by 11 intervals. The concentration of the activating molecule is measured at the activating site and the inhibitory molecule at the inhibitory site.

By analyzing curves similar to those in figure 3 the frequency of oscillations was calculated. Figure 4 illustrates the sudden appearance of high frequency oscillations when the sites are near each other and the decreasing frequency as they are moved apart until oscillations disappear.

Figure 4 - Frequency of Activation Molecule at Activation Site vs Separation Distance – The symmetrical appearance is a result of the activation site at first moving away from the inhibitory site, making its way around the ring, and then getting closer to the inhibitory site along the opposite side. The graph of frequency of the inhibitory molecule at the inhibitory site vs separation distance was equivalent.

The observation that the steady state becomes unstable for a range of separations has been reported previously by Shymko and Glass (1974).

Another way of looking at how oscillations depend on separation distance is by plotting the value of a curve whenever its derivative is zero vs the separation distance. This will also tell us some information on how the shape of the oscillations change with separation distance.

Figure 5 - Amplitude of Oscillations of Activation Molecule at Activation Site vs Separation Distance – The amplitude of oscillations was generated by noting the value of the curve whenever its derivative was zero and plotting this for each separation distance. So that separation distances that have only one value plotted are at steady state. Whereas separation distances with two values plotted are oscillating. The larger of the two values corresponds to the peak and the smaller the trough.

Notice the range of distances that have two values plotted correspond to the range of distances that have a nonzero frequency in figure 4.

One concern is that it is unlikely the chosen parameters are representative of the synthetic genetic network that we are constructing. The major parameter that we currently have no control over is the Hill coefficient n. We believe that if the Hill function is steep enough then the circuit will act as an on/off switch and oscillations are possible. Increasing n will only increase the steepness of the Hill function so we are safe in assuming that n > 8 will still lead to oscillations. However what is the smallest n that would still allow oscillatory behaviour? We calculated frequency for the range of separation distances with the following Hill coefficients n=2,3 and 4.

Figure 6 - Frequency of Activation Molecule at Activation Site vs Separation Distance (n=4) - Oscillations still present for a similar range of separation distances as n=8 but with lower frequencies.

Figure 7 - Frequency of Activation Molecule at Activation Site vs Separation Distance (n=3) - No oscillations observed.

In order to observe oscillations we require a Hill coefficient of at least n=4.

We also looked into the dependence on the diffusion rate, D, and the presence of oscillations. As can be imagined, if the diffusion rate is too slow or too fast oscillations will not occur. This was confirmed with preliminary simulations where we found a particular separation distance that loses oscillatory behaviour if the diffusion rate is outside a critical range. We endeavoured to form a bifurcation diagram showing the ranges of separation distance and diffusion rate where oscillations occur, however this task is computationally intensive and we have not completed it as of this writing. We hope to have this result in time for the Jamboree.

Two Oscillators

We next looked at a system consisting of two oscillators, where each consists of an activation and inhibitory site.

Figure 8 – Two Oscillators – Each oscillator consists of one activation and one inhibitory site. There are two arrangements of the system: BA AB, termed the symmetrical system, and BA BA, the unsymmetrical system. In this document we only discuss the symmetrical system. Oscillator 2 travels around the ring.

The distance between the two oscillators was varied while the distance between two sites within an oscillator was held fixed at 5 intervals. This value was chosen for demonstration purposes, however the dynamics to be described have been observed at various separation distances.

We first looked at the change in frequency as the oscillators are moved apart.

Figure 9 - Frequency of Activation Molecule at Activation Site of Oscillator 2 vs Separation Distance – There are eight separation distances that resulted in a frequency of -1: 14,15,16,17,474,475,476, and 477. This value is used as a marker to denote a curve who's frequency did not converge to a single value. As will be seen later, these distances resulted in complex dynamics.

The first striking feature of this graph is the appearance of separation distances that have negative frequencies. The value of -1 is assigned to a separation distance that results in a complex oscillation whose frequency cannot be resolved using a simple method, which assumes simple periodic motion. Before we look into these curves, let's look over the entire graph. When the two oscillators are very close together, their frequency is always less than an isolated oscillator (look at figure 4 for two sites with a separation distance of 5 intervals). Interestingly, between the separation distances of 1 and 10 there appears to be a local minimum frequency. Meaning as the oscillators are moved apart, their frequency initially decreases and then begins increasing. This trend of increasing frequency continues until the numerical method calculating frequency breaks down. To get a better idea of why the frequency calculating method is not converging let's look at the amplitude of oscillations vs separation distance:

Figure 10 - Amplitude of Oscillations of Activation Molecule at Activation Site of Oscillator 2 vs Separation Distance - Notice that between the distances 14-18 there appears to be more than 2 amplitude values per distance.

We can now begin to understand why the frequency calculating method failed. This method assumes simple periodic motion (similar to a sine or cosine wave), however the curves for which a frequency could not be calculated have more than 2 values where their derivative equals zero. This hints that they are not simple periodic curves. The following figure illustrates the concentration of the activation molecule at the activation site of oscillator 2 for the first distance that results in an oscillation whose frequency cannot be computed.

Figure 11 - Concentration of activation molecule at the activation site of oscillator 2 for a separation distance of 13 - There appears to be a constant "resetting" of the oscillations. As if a threshold is constantly being crossed that induces a reset of the oscillatory behaviour.

It is immediately obvious why the frequency measuring numerical method failed, the oscillations are far from simple. Before investigating this curve further, let's take a look at the other separation distances that resulted in complex oscillations.

Figure 12 - Concentration of activation molecule at the activation site of oscillator 2 for a separation distance of 14

Figure 13 - Concentration of activation molecule at the activation site of oscillator 2 for a separation distance of 15

Figure 14 - Concentration of activation molecule at the activation site of oscillator 2 for a separation distance of 16

The separation distances 474,475,476, and 477 resulted in similar dynamics and are not shown for brevity.

We first wanted to know whether these were truly the only separation distances that resulted in this type of dynamics. We recalculated the frequency of oscillations at each separation distance after simulating for 10,000 time units. In this way we hoped to overcome any transient dynamics and observe the true interaction. This resulted in the same set of separation distances having complex oscillations while everywhere else had regular simple oscillations. Furthermore, we recalculated the frequency of oscillations at each separation distance by using the final time point of the immediately preceding separation distance as the initial condition. We thought this would allow a quicker convergence to the stable limit cycle. However this also resulted in the same findings as before.

In order to further classify these dynamics as either aperiodic, quasiperiodic, or chaotic we must look further than the time evolution graphs. The phase plots begin to describe what type of oscillations are present. We plot the concentration of the activation molecule at the activation site of oscillator 1 vs the activation molecule at the activation site of oscillator 2.

Figure 15 - Activation molecule at the activation site of oscillator 1 vs the activation molecule at the activation site of oscillator 2 for distance of 15 - This was compiled from a simulation lasting 100 time units.

Figure 16 - Activation molecule at the activation site of oscillator 1 vs the activation molecule at the activation site of oscillator 2 for distance of 17 - This was compiled from a simulation lasting 100 time units.

Figure 17 - Activation molecule at the activation site of oscillator 1 vs the activation molecule at the activation site of oscillator 2 for distance of 17 - This was compiled from a simulation lasting 0.156 time units. Notice the colours correspond to numerical simulation time not real time.

Comparing figure 11 to figure 12 lets us believe that the oscillations in figure 11 are periodic while those from figure 12 are quasiperiodic. This is further illustrated in figure 13 where we can see that on each cycle there is a small shift in the oscillations, which is characteristic of quaisperiodic curves. In order to further investigate these dynamics we look at the Poincare maps of each curve.

Figure 18 - Poincare map for separation distance of 15 intervals - Since the number of dots do not increase with more iterations we are convinced that this is a periodic oscillation.

Figure 19 - Poincare map for separation distance of 17 intervals - The number of dots increase with each iteration, meaning this is truly a quasiperiodic motion.

Since the Poincare map of the oscillations present at a separation distance of 17 intervals leads to a filled in closed shape we are fairly confident that we have found quasiperiodic motion. Although the curve seen at a separation distance of 15 also looks quasiperiodic, the Poincare map does not form a filled in closed shape we can conclude that this separation distance leads to complex oscillations but not quasiperiodic.

The final analysis that we will perform on these simulations is to look at the phase lag between the two oscillators. We again move the second oscillator around the ring while measuring the phase lag between the activation molecule at the activation site of the first oscillator and the activation molecule at the activation site of the second oscillator.

Figure 20 - Phase lag between the activation molecule at the activation site of the first oscillator and the activation molecule at the activation site of the second oscillator - The phase lag is calculated by finding the lag in real time and then dividing by the period of oscillations. Negative phase lag values correspond to those separation distances where a frequency could not be calculated.

A striking result of the phase lag calculation is that for some range of distances past the quasiperiodic forming distances, all the oscillations are exactly 180 degrees out of phase! A phase plot gives us more insight into these dynamics.

Figure 21 - Activation molecule at the activation site of oscillator 1 vs the activation molecule at the activation site of oscillator 2 for distance of 20 - The perfect reflection across the line y=x means the two curves are always 180 degrees out of phase.

It is currently not well understood why a range of separation distances past the quasiperiodic curves remain 180 degrees out of phase.

Two Dimensional Modeling

We were able to further observe oscillatory behaviour in two dimensions, unfortunately the .avi files we have are too large to upload onto the wiki. This is all we have right now to describe our efforts into two dimensional modeling. We'll hopefully be describing them in detail during our presentation at the Jamboree.

Appendix A - Parameters

The standard set of parameters used to observe oscillations were taken directly from Shymko and Glass (1974).

γ = 2, D = 2, λ = 54, θ = 1, b = 0, N = 8.

Both strains were assumed to have identical parameters in order to simplify the model as well as explicitly observe the dependence of dynamics on separation distance.

• Recent changes
• What links here
• Related changes
• Special pages
• My preferences

• Printable version
• Permanent link
• Privacy policy
• Disclaimers