ENT 420 Biological System Modeling Phase Plane Diagrams

ENT 420 Biological System Modeling Phase Plane Diagrams Lecturer Engr. Mohd Yusof bin Baharuddin MBiomed. Eng (Melbourne) BBiomed. Eng (Malaya)

ENT 420 Biological System Modeling Objectives • • Motivation Graphical methods for exploring systems of ODEs Phase plane diagrams Example/tutorial: Chemostats Numerical software: xppaut Limit cycles, oscillators Introduce the FH-N class of models The Hopf Bifurcation

ENT 420 Biological System Modeling The story so far…. Model type Difference equations Systems of ODEs General form Steady state condition Stability condition Graphical exploration Cobweb graph Phase plane diagrams

ENT 420 Biological System Modeling Graphical interpretation of ODEs Consider When y=0 and y=1, y 1 0 t

ENT 420 Biological System Modeling Graphical interpretation of ODEs When 0<y<1, When 1<y,

ENT 420 Biological System Modeling Graphical interpretation of ODEs y 1 0 t

ENT 420 Biological System Modeling What about systems of ODEs? Too hard! y We could produce a 3 D plot i. e. with axis (x, y, t) y x t t x

ENT 420 Biological System Modeling Graphical interpretation of ODEs Instead let’s look at how the slope and direction of the vector (dx/dt, dy/dt) varies within a plane (x-y 2 D plot) e. g. y y x A phase plane graph x Slope of vector given by

ENT 420 Biological System Modeling Phase plane graph • Example: Recall the Lotka-Volterra predator-prey model prey Oscillatory neutral stable steady state predator 1

ENT 420 Biological System Modeling Phase plane graph • Instead of graphing the behaviour of both predator and prey populations as a function of time, we instead can look at the change in prey populations versus change in predator populations on a phase plane plot. • That is…. .

ENT 420 Biological System Modeling Phase plane graph • Rearrange to and integrate rearrange

ENT 420 Biological System Modeling Phase plane graph prey v predator 1 c 4 1 0 c 3 c 2 u 1 Neutrally stable Can’t always integrate to get exact solution like in this example.

ENT 420 Biological System Modeling A more systematic approach to phase plane diagrams Draw the nullclines: Let’s look at a previous example y X - nullclines y-nullclines y=1 steady states x=1 x-nullclines x y - nullclines These define the nullclines Steady states (singularities) occur at the intersection between x nullclines and y nullclines

ENT 420 Biological System Modeling A more systematic approach Assign direction along the nullclines according to vector: Note one component of the vector will be zero on a nullcline. Arrows on y nullclines parallel to x axis Arrows on x nullclines parallel to y axis

ENT 420 Biological System Modeling A more systematic approach We can start to see the flow patterns already Nullclines separate the phase plane plot into regions in which the flow is in the same general direction Dramatic local changes in the flow pattern can really only occur in the vicinity of steady states Begin to see global behaviour i. e. behaviour not local to steady states

ENT 420 Biological System Modeling Flow patterns around (linear) steady states in phase plane diagrams Image from Fall et al. “Computational Cell Biology” (2002).

ENT 420 Biological System Modeling Murray (2002)

ENT 420 Biological System Modeling Gorillas again Identify and draw nullclines H H=k/m 0 N=0 H=r. N/j These define the nullclines Steady states occur at the intersection between x N nullclines and y nullclines steady state

ENT 420 Biological System Modeling Gorillas again Assign flow direction arrows along nullclines When N=0 H When H=k/m 0 N=0 H=r. N/j When H=r. N/j N steady state

ENT 420 Biological System Modeling Gorillas again What kind of steady state do we have? Can’t tell just from phase plane plot! H H=k/m 0 N=0 H=r. N/j N steady state

ENT 420 Biological System Modeling k=1. 1 m 0=0. 1 at t=0 N=1000 and H=5 r=0. 01 j=0. 1

ENT 420 Biological System Modeling Tutorial Example: Chemostat Bacterial growth in a Chemostat A chemostat is a device for harvesting bacteria. Stock nutrient of concentration c 0 enters the bacterial culture chamber with a flowrate F. For mass conservation there is an equal flow rate F out of the culture chamber. V is the volume of the culture chamber and N is the number density (number per volume) of bacteria in the culture chamber F c 0 N c V F

ENT 420 Biological System Modeling Bacterial growth in a Chemostat The equations for the rate of change of the nutrient concentration and bacteria number density are as follows: F c 0 c N V F k(c) is the concentration dependent reproduction rate of the bacteria α is units of nutrients consumed in producing one unit of population growth Note F has units of volume/time, FN/V has units number density of bacteria per time.

ENT 420 Biological System Modeling Example: Chemostat Bacterial growth in a Chemostat The rate of growth of bacteria increases with nutrient availability only up to some limiting value. (The poor little individual bacterium can only consume nutrient and reproduce at some limited rate. ) k(c) kmax Michaelis-Menton kinetics (more next week or the week after) 0. 5 kmax c 1/2 c

ENT 420 Biological System Modeling Bacterial growth in a Chemostat Non-linear ODE: We can non-dimensionalise these equations to make them look neater and to reduce the number of parameters. Let t=t*V/F, c=c*c 1/2 and N=N*c 1/2 F/(αVkmax), where N*, C* and t* denote dimensionless bacteria density, nutrient concentration and time. We can then write the ODEs as: where α 1=Vkmax/F and α 2=c 0/c 1/2

ENT 420 Biological System Modeling Example: Chemostat For simplicity of mathematics let’s drop the superscript * and let’s assume that α 1=2 and α 2=3

ENT 420 Biological System Modeling Tutorial: Phase plane diagram of chemostat model c Find nullclines: N nullclines -> c nullclines -> 0, 0 N

ENT 420 Biological System Modeling Example: Chemostat Bacterial growth in a Chemostat Non-linear ODE: Steady states:

ENT 420 Biological System Modeling Example: Chemostat Bacterial growth in a Chemostat Stability of steady states: Linearise our equation to obtain the A matrix: Recall where

ENT 420 Biological System Modeling Example: Chemostat Bacterial growth in a Chemostat Stability of steady states: Linearise our equation to obtain the A matrix: Find eigenvalues at each steady state by

ENT 420 Biological System Modeling Example: Chemostat Stability of steady states: When

ENT 420 Biological System Modeling Example: Chemostat When

ENT 420 Biological System Modeling Example: Chemostat

ENT 420 Biological System Modeling Software • XPPAUT (freeware) (http: //www. math. pitt. edu/~bard/xpp. html) • a tool for solving – – – differential equations, (up to 590 differential equations) difference equations, delay equations, functional equations, boundary value problems, and stochastic equations. • Also creates phase plane diagrams, cobweb graphs does stability and bifurcation analysis etc (unlike Matlab)

ENT 420 Biological System Modeling Xppaut Program • Write program in a text editor • The following is program for a linear system of 2 ODEs. # equations dx/dt=a*x+b*y dy/dt=c*x+d*y # parameters par a=0, b=1, c=-1, d=0 # intial conditions init x=1, y=0 # that's all folks done Minimum program dx/dt=a*x+b*y dy/dt=c*x+d*y par a, b, c, d done

ENT 420 Biological System Modeling Xppaut Program • Write program in a text editor • The following is program for the Lotka-Volterra system prey du/dt=u*(1 -v) dv/dt=a*v*(u-1) par a=1 done predator 1 v 1 0 1 u

ENT 420 Biological System Modeling Xppaut Program • Write program in a text editor • The following is program for the system considered at the start of this lecture du/dt=u*v-v dv/dt=v*u-u done

ENT 420 Biological System Modeling Xppaut Program • Write program in a text editor • The following is program for the chemostat c d. N/dt=alpha 1*(c/(1+c))*N-N dc/dt=-(c/(1+c))*N-c+alpha 2 par alpha 1=5, alpha 2=0. 5 @ xp=N, yp=C, xlo=-. 25, xhi=3, ylo=-. 1, yhi=1, total=100 @maxstor=10000 done N

ENT 420 Biological System Modeling Xppaut Program • Write program in a text editor • The following is program for the cancer model in lecture 5 dc/dt=delta*c*(w*c/(1+I)-1) d. I/dt=sig+b*(gamma*c-I) par delta=1, W=1, sig=0, b=1, gamma=0. 2 @ xp=c, yp=I, xlo=0, xhi=3, ylo=0, yhi=3, total=100 done

ENT 420 Biological System Modeling At critical points S, R, Q the vector V is zero, and corresponds to stationary points Closed paths like C correspond to periodic solutions Simmons, G. (1991) “Differential equations with applications and historical notes”

ENT 420 Biological System Modeling Poincaré-Bendixson theorem If the phase plane contains a domain, enclosed by a boundary B on which the vector always points into the domain, and the domain contains a singular point P which is an unstable spiral or node then any phase trajectory cannot tend towards the singularity or leave the domain. Murray (2002) The Poincaré-Bendixson theorem says that as trajectory will tend to a limit cycle the

ENT 420 Biological System Modeling Oscillators, excitability and Fitz. Hugh-Nagumo Models • Fitz. Hugh and Nagumo independently proposed a simplified model of the Hodgkin-Huxley equations that are used to describe the variation in sodium and potassium ions across a cell membrane of a nerve cell (neuron). • Now, FH-N models are those with one linear nullcline for the “slow variable” and a cubic nullcline with an inverted “N” shape for the “fast variable”

ENT 420 Biological System Modeling Oscillators, excitability and Fitz. Hugh-Nagumo Models Example: x represents the excitation variable “fast variable” y represents the recovery variable “slow variable” Parameters a, b and c in the model are assumed to be positive with the following restrictions z is a stimulus intensity

ENT 420 Biological System Modeling Phase plane diagram of FH-N model y z=0 x Two Nullclines: y z=-0. 4 x Slight shift in curve

ENT 420 Biological System Modeling Oscillators, excitability and Fitz. Hugh-Nagumo Models Stability? For the trace to be negative and the determinant positive i. e. a stable steady state or

ENT 420 Biological System Modeling Oscillators, excitability and Fitz. Hugh-Nagumo Models That is the steady state is not stable if it falls within the region Imaginary number Since the requirement for a real value steady state not being stable is just Note: The stability condition of the steady state does not change with z, but the steady state itself does.

ENT 420 Biological System Modeling Xppaut Program • Write program in a text editor • The following is program for the FH-N model du/dt=c*(v+u-(1/3)*u^3+zed) dv/dt=-(1/c)*(u-a+b*v) par a=0. 7, b=0. 8, c=3, zed @ xp=u, yp=v, xlo=-3, xhi=3, ylo=-3, yhi=3, total=100 @maxstor=10000 done

ENT 420 Biological System Modeling

ENT 420 Biological System Modeling The Hopf Bifurcation • Hopf Bifurcation Theorem predicts the appearance of a limits cycle about any steady state that undergoes a transition from a stable to an unstable focus (spiral) as some parameter is varied.

ENT 420 Biological System Modeling The Hopf Bifurcation Consider A steady state occurs at At this steady state

ENT 420 Biological System Modeling The Hopf Bifurcation du/dt=v dv/dt=-v^3+r*v-u par r @ xp=u, yp=v, xlo=-3, xhi=3, ylo=-3, yhi=3, total=100 @maxstor=10000 done

ENT 420 Biological System Modeling Summary ØPhase plane diagrams, nullclines, shapes of flow patterns near various linear steady states ØExamples: Chemostat, gorillas ØXppaut ØLimit cycles ØExcitability ØFH-N models

ENT 420 Biological System Modeling Somewhere people still eating once a day…. ~Yusof~