Team:Aberdeen Scotland/modeling/pde

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.

Where 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.

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 quorum sensing, nonlinear diffusion and cell-growth & multiplication. The following documents the theoretical knowledge necessary to create a numerical scheme to solve the system - which may be undertaken on a "Next steps" basis.

Comment is made on the extension of the Minimal Model to incorporate quorum sensing and cell-growth & multiplication - so as to be fully representative of the PicoPlumber system - in the final subsection of this page.

Minimal Model

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

where is the chemotactic sensitivity, the bacterial density, the concentration of the chemoattractant and , their respective time derivatives.

In Ephsteyn & 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.

To solve the Minimal Model we require a Mesh upon which to discretise the system and a reformulation of the Convection-Reaction-Diffusion to render it soluble by Discontinuous Galerkin methods.

The following documents the theoretical considerations necessary to obtain a robust solution. Ephsteyn & Kurganov^[4] must be fully credited for their description of the numerical scheme outlined below - which is expanded upon concomitantly with the considerations of Hillen & Painter^[1]- in a manner so as to provide clarity and insight to the reader.

Mesh Generation

The method begins by defining a Sobolev space over the domain 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

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

And for those w in the Sobolev space such that for w evaluated at an arbitary element ; w is in the space of polynomials of degree over , styled

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

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 , we obtain the variation

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]

where the boldface 3-vectors are , , and .

Ephsteyn & Kurganov^[4] however, report that such a re-formulation is prone to severe instabilities in the transition from a Hyperbolic region to an Elliptic one - upon implementation of a high-order Numerical Method. To implement the Discontinuous Galerkin method we therefore rewrite the system as

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

once again where the boldface 3-vectors in this instance are where , are the fluxes and is the reaction term. It should similarly be noted that k=1 is a constant and the convective part of the system is requisitely Hyperbolic.

This renders the Minimal Model soluble by a Discontinuous Galerkin numerical scheme as proposed by Ephsteyn & Kurganov^[4].

Together, the Mesh and the Convection-Reaction-Diffusion equation form the complete theoretical basis upon which the program is to be developed.

Analytic Form of Solution

In the nomenclature defined above, a solution to the Keller-Segel System must take the form

for some solution in the direct product space of the respective nondegenerate quasi-linear rectangular subdivisions , satisfying the weak formulation of the Minimal Model

where sums of the integrals are over and the set of interior edges , which in turn are split up into horizontal & vertical edges - and respectively - of the rectangular boundary .

Feasibility of developing Numerical Method

As seen in the precending section. even the simplest of analytic solutions to the Keller-Segel System using the Discontinuous Galerkin method is distinctly non-trivial - and would require a high degree of technical proficiency and time to implement. This was outwith the scope of our project, but were it to be attempted, would prove a robust and highly accurate means of collimating the individual PicoPlumber modelling modules.

Moreover, comment must be made on extensions to the basic Minimal Model "chassis" to incorporate effects such as quorum sensing and cell-growth & multiplication. These augmentations are discussed thoroughly in Hillen & Painter^[1] but would require a re-derivation of some of the results given on this page - assuming compatibility with the Discontinuous Galerkin methods, that is.

Back to Hill Input Functions

Continue to Chemotaxis

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