Difference between revisions of "Team:Freiburg/Results/Modeling"

It is an inhomogeneous parabolic partial differential equation (PDE). The diffusion constant D depends on the media and materials involved in the system, the inhomogenity describes sources of the system. The velocity vector v is assumed to be a possible sedimentation velocity (A4) and therefore a constant. The initial and boundary conditions are the following: The PDMS slide as well as the iRIf slide are closed surfaces while the left and right edges of the system are not. On the PDMS slide proteins are produced by a logistic function (A5). The solution to this system is only numerically solvable. Therefore it has to be prepared mathematically by discretizing the whole setup. Finally we used the Crank-Nicholson approach on a finite-differences-method to solve this system.

3.3 Numerical approach

Let us start by describing the crosssection in terms of mathematics: It is a subset M of the set of twodimensional real vectors. Of this subset we define two more subsets S and B describing the DNA spot and the glass slide respectively.

Two time &\delta; is the spot width, its thickness is ε. The thickness of the binding area is ξ yet we choose them to be the same (A6). Next we have to simplify our PDE: It consists of three terms on the right side which can be (in our case) shortened to the following components:

The next step is the change from a continuous to a discrete system: First we describe our concentration function as one depending on three integers (places x and z, time t):

Also, instead of handling every component individually we define a concentration vector containing all concentration values of M at a static time step:

We do the same for our source by defining a production vector at the time t:

Why is this done? Because we can split our grid into a connected line of rows. The first row of the grid then would be given by the first components of the vector, followed by the second row and so on.

Next we discretize the derivations of the PDE by exchanging the differentials by finite differences. In this case differentials can be seen as the rate of change in one (x-,z- or) t-step (secants so to say). For a first and second derivation this would be written as:

Using the integer multiplicated step values for our PDE components this results in:

Logically this method is called "finite differences method". For conveniences sake two parameters λ and μ are defined and used in the upcoming simplifications:

There are different approaches with the finite differences method to approximately solve the PDE. We decided to use an ansatz by Crank and Nicholson which has the following form: A functions change in one time step is determined by the average of momentary and following time step of the right hand side of the PDE. Or, in mathematical terms:

Combining the PDE parts and splitting the concentration values into following and momentary time steps the problem can be solved by solving a linear system of coupled equations in place components of the concentration:

This is an iterative method: The values of the upcoming time step are calculated by the ones from the former time step. The equations cannot be solved independently from one another. In a grid of NW rows and NL columns this a system of NW times NL equations that have to be solved every step.

What can be seen in the scheme above is the underlying geometrical ansatz: Diffusion and convection are done only with direct grid neighbours (two vertical and two horizontal ones). Logically this holds true only for non-border grid points. For the systems borders our boundary conditions have to be considered which manifest in the following conditions:

Inserting those in the main Crank-Nicholson-scheme leads to the border descriptions:

Combining two border conditions leads to corner conditions which have to be considered too. They follow the same calculation yielding:

Being a linear system of equations this problem can be replaced by a short form using matrix-vector multiplication (hence the definition of the vectors from before):

Where the two matrices A and B consist of the factors accompanying the concentration values. They are rather large and are sketched below:

This concludes the theoretical part of the diffusion. The next step was to implement this system using python 3.4.2 including numpy and matplotlib. Source codes can be obtained here:

The first file produces two output files. One of them contains the protein distribution on the glass slide at the finishing time. The other one consists of the total bound protein at the surface of the slide. The second file returns pictures of the ongoing diffusion giving an insight of the diffusion process in the crosssection system.

3.3 Assumptions

A1: The spots are equally distributed and distance to the borders of the microfluidic chamber is at least one spot-to-spot distance.

A2: The spots are perfect circles and the binding site at the iRIf glass is homohegeneous -> Around one spot cylindrical symmetry is given.

A3: The area of produced proteins is as thick as the one of bound ones -> χ = ε

A4: Movement by virtue of gravitation is considered to already be in a stationary status (constant movement).

A5: Protein production is described by a logistic function with a minimal value close to zero and maximum value of 100 (percent).

Results

1. Transcription and Translation

We wanted to test our model's ability to predict the concentrations of the cell-free expression biochemical network entities and to adapt to different conditions. Given this properties, our model has the potential to reveal limiting resources and bottleneck reactions and to test optimized parameter settings for the resulting effects on the biochemical network behavior.

We run the model using the default starting concentrations and parameter values as given in the following, tracking the concentration time course of free RNA polymerase, free 30s and 50s ribosomal subunits, free 70s ribosome, free amino acids, free tRNA, free initiation factors, free elongation factors, free nucleotides, free mRNA and free protein, to examine the network's predicted behavior. Afterwards, we tracked the protein and mRNA concentrations in three more simulations, changing one entity starting concentration each, to examine the network's predicted reaction and compare it to the reaction expected.

Default Parameter Values

Network Behavior Under Default Parameter Values

Network Reaction To Changed Parameter Values

Summary and Model Optimizations

Our model describing cell-free expression, including 65 entites and 44 parameters, showed the ability of both predicting network behavior and network reaction to changed conditions. Although the model's predicted output does not fit to real cell-free expression yet, there is a huge capability of making the model more specific. Using the set of assumptions defined, the model can be extended step by step, bringing it closer to the cell-free expression biochemical network.

3. Diffusion results

For once we made an animation of the diffusion using exemplary values. You can see a video of the diffusion process down below:

Figure: ...

Figure: ...

Also we checked the protein binding for different spot widths: 20%, 40%, 60% and 80% of the total crosssection width. While the binding distribution returned gaussian profiles the total binding like the production