Team:Aberdeen Scotland/modeling/pde

From 2009.igem.org

Revision as of 19:33, 21 October 2009 by Takemusu (Talk | contribs)

University of Aberdeen iGEM 2009

Contents

Introduction

It was realised quite soon into the project that upon modelling the Inner Dynamics of the system, the individual modules could be collimated in to a single Partial Differential Equation – hereafter referred to as PDE – which would accurately describe the evolution of the system in space and time.

A PDE model was an excellent “next step” for our modelling effort – and although a physically meaningful model could not be generated; the following describes a concise modelling methodology which would merit implementation outside of iGEM – or indeed, in future iGEM competitions.

The Keller-Segel Chemotaxis system is a set of two coupled Nonlinear PDE’s; ostensibly Parabolic except for a mixed Hyperbolic-Elliptic term in the cross-diffusive flux. The terms “Parabolic”, “Hyperbolic” and “Elliptic” refer to the number of boundary conditions the system must satisfy – the dependence on which defines (or “sets”) the behaviour of the system.


PDE1.gif

Where Var34.gif is the cross-diffusive flux

The Keller-Segel system was perfect for modelling our PicoPlumber modules, inasmuch as Quorum Sensing could be effectively integrated into the formulation[1]. This allowed the advantage of analysing the dynamical evolution of the system in space & time - and in regard to the Internal Dynamics.

The following is an account of our researches into utilising these powerful Mathematical tools to describe our system.

One-Dimensional Diffusion Equation

As PDE’s are a challenging and actively-researched topic in Mathematics[2-3]; our modelling effort began by solving a much simpler problem by way of introduction to the field.

The one-dimensional Heat Equation is the classic introduction to Partial Differential Equations. It was modelled by utilising an Implicit Finite Difference method of our own composition. This allowed us to attain an aptitude for tackling the more difficult Keller-Segel system to come.

The Heat Equation and the Diffusion Equation of our chemo-attractant shared the same mathematical formulation – differing in fact, by the choice of the constant k – and so this introduction still retained a Biological pertinence.


PDE implicit.gif

Keller-Segel System for Chemotaxis

As mentioned above, the Keller-Segel System for Chemotaxis is the set of coupled Nonlinear PDE's which describe chemotactic motion. The mixed Hyperbolic-Elliptic term endows the system with rich behaviour – but also renders it insoluble by elementary Numerical Methods[4].

In searching for a more advanced numerical method, Discontinuous Galerkin Finite Elements were found in Systems Biology literature[4]. to have successfully modelled the system and identified parametric instabilities.

In their award-winning paper[5]. , Hillen and Painter describe augmentations to the standard Keller-Segel System, hereafter referred to as the "Minimal Model". These augmentations account for effects such as density- and signal-dependent sensitivity, nonlinear diffusion and cell-growth & multiplication. The following documents the theoretical knowledge necessary to create a numerical algorithm to solve the system - which may be undertaken on a "Next steps" basis.


Minimal Model

The Chemotactic Minimal Model is the simplest "chassis" for describing a bacterial population chemotaxing towards a chemoattractant. It is of the form

PDE1.gif

where Var3.gif is the chemotactic sensitivity, Var5.gif the bacterial density, Var4.gif the concentration of the chemoattractant and Var1.gif, Var2.gif their respective time derivatives.

In Ephsteyn and Kurganov[4], the aforementioned Discontinuous Galerkin Finite Element Method is applied to successfully solve the system. The following will provide an abridged summary of the method, covering all elements necessary to implement it.

Theoretical Considerations

THe Finite Element Method is a numerical technique for the solution of PDE's by discretising the domain over which the PDE is solved into finite "elements". The discretised domain is then termed the mesh. Elements are interpolated by polynomial basis functions which serve to approximate the equation over the element without allowing errors in intermediate calculations to accumulate and render the final output meaningless.

Discontinuous Galerkin methods are an augmented family of Finite Element Methods developed by the Russian Engineer, Boris Galerkin. They are superior to other finite element methods for several reaons, including - but not restricted to - local, element-wise mass conservation and flexibility to choose non- and high order polynomial bases. Element-wise mass conservation is important, and naturally the bacteria must be able to "exist" computationally through all space & time.

Mesh Generation

The method begins by defining a Sobolev space Var14.gif over the domain Var31.gif upon which the mesh is specified. A Sobolev space is a Linear Space of polynomial functions - strictly a Banach space of square-integrable functions - endowed with the Sobolev s-norm

PDE16.gif
,
PDE17.gif


For s a labelling-index and Var16.gif an element of the discrete space of discontinuous piecewise polynomials of degree Var32.gif as a function of the nondegenerate quasi-uniform rectangular subdivisions Var13.gif

PDE15.gif

And for those w in the Sobolev space such that for w evaluated at an arbitary element Var33.gif; w is in the space of polynomials of degree Var32.gif over Var16.gif, styled Var15.gif

While this may seem mathematically "verbose", the mesh itself and the non-degeneracy of the elements therein can be seen in the diagram of Var35.gif below - whilst the Linear space constructions remain abstract

Fea mesh2.gif

This rigourously defines the mesh upon with the Discontinuous Galerkin method solves the Keller-Segel system - which is the key component of the algorithm.

Convection-Reaction-Diffusion Reformulation

Algebraically manipulating the Minimal Model by defining the vector Gradc.gif, we obtain the variation

PDE7.gif

which is in a form soluble by Discontinous Galerkin methods. This form of the Minimal Model can in turn be considered as a Convection-Reaction-Diffusion system[4]

PDE9.gif

where the boldface 3-vectors are U matrix.gif, Fu matrix.gif, Gu matrix.gif and Ru matrix.gif. Ephsteyn however, reports that such a re-formulation is prone to severe instabilities in the transition from a Hyperbolic region to an Elliptic - upon implementing a high-order Numerical Method. To implement the Discontinuous Galerkin method we therefore rewrite the system as

Chem system dg.gif

subject to the boundary conditions PDE13.gif. This may also be considered as a Convection-Reaction-Diffusion system of the form

PDE14.gif

once again where the boldface characters Dg sys q.gif, Dg sys fq.gif, Dg sys gq.gif and Dg sys rq.gif are the fluxes, with k=1 a constant - where the convective part of the system is Hyperbolic.

References

[1] Hillen T, Painter KJ, "A user's guide to PDE models for Chemotaxis" J Math Biol. 2009 Jan;58(1-2):183-217. Epub 2008 Jul 15.

[2] Oxford Centre for Nonlinear PDE. http://www.maths.ox.ac.uk/groups/oxpde/

[3] The Maxwell Institute: Centre for Analysis and Nonlinear PDE's. http://www.maxwell.ac.uk/PDEs/

[4] Ephsteyn Y, Kurganov A, New Interior Penalty Discontinuous Galerkin Methods for the Keller–Segel Chemotaxis Model. J. Numer. Anal. Volume 47, Issue 1, pp. 386-408 (2008)

[5] Faculty of 1000 Biology. "A user's guide to PDE models for Chemotaxis" Importance recognition. http://f1000biology.com/article/id/1126950