### Model-based Design in Synthetic Biology - Mathematics

```Mathematical Modelling for
Synthetic Biology
aGEM Workshop on Mathematical Modelling
July 22, 2012, Lethbridge, AB
Brian Ingalls
Department of Applied Mathematics
University of Waterloo
Waterloo, ON
Workshop Outline
• Introduction to mathematical
modelling of biochemical reaction
networks
• Modelling of gene regulatory
networks
• Lab I: simulation of kinetic models
• Tools for model analysis
• Lab II: model-based design of
gene regulatory networks
Models in Science and
Engineering
• Models are abstractions of reality
Models in Science and
Engineering
• Models are abstractions of reality
• Models can be physical
Ball-and-stick model of
molecular structure
http://mariovalle.name/ChemViz/representations/index.html
Mouse model of obesity
http://srxawordonhealth.com/2012/04/
Models in Science and
Engineering
• Models are abstractions of reality
• Models can be physical or conceptual
Interaction diagram model
of G-protein signalling
http://www.nature.com/scitable/topicpage/gpcr-14047471
Kinetic model of bacterial
chemotaxis signalling pathway
Models in Science and
Engineering
• Models are abstractions of reality
• Models can be physical or conceptual
• Mathematical models are mechanistic
(based on physico-chemical laws) and
predictive (allow inferences beyond the
data used for their construction)
How are mathematical models
used in molecular biology?
• Models summarize data
• Models allow of falsification of hypotheses
• Models allow exploration of system
behaviour (in silico experiments)
• Model-based design allows easy exploration
of design space
Model Construction
Modelling Chemical Reaction
Networks
Chemical reaction:
Rate constant:
Law of mass action:
Using derivatives to describe
rates of change
Ex: decay reaction:
differential
equation
model
rate of change of [A] at time t
rate of reaction at time t
Solution of the differential
equation model
Model simulation = in silico experiment
Model Simulation
Numerical simulation of
differential equation models
Approximate derivative by a difference quotient
Rearrange to yield an update rule:
Repeated application of the update rule (starting
from known initial concentration):
Implemented in MATLAB, XPPAUT, Copasi,
Mathematica, Maple, …
Tutorials in notes for
XPPAUT (freeware, simulation and analysis of differential equation models)
MATLAB (licensed, general-use computational software)
Example network model
networ
k
model
simulation
Model Analysis
Separation of time-scales
Every model is formulated around a specific time-scale
 Processes occurring on a slower time-scale are treated as frozen in time
 Processes occurring on a faster time-scale are presumed to occur
instantaneously
Treating rapid processes
Rapid equilibrium approximation: presume
at all times.
Quasi-steady-state approximation: presume [A] is in steady-state with
respect to [B] at all times.
Sensitivity analysis
Measure sensitivity of steady-state species concentrations to changes in
model parameters
concentration
s (p=p1)
}ds
s (p=p2)
time
Applications of Sensitivity
analysis
Identification of optimal
drug targets (steps with
high sensitivities)
Bakker et al. 1999, ‘What controls
glycolysis in bloodstream form
Trypanosoma brucei? JBC 274.
Applications of Sensitivity
analysis
Interpretation of regulation schemes:
role of negative feedback
Biochemical Kinetics
Saturation: Michaelis-Menten
kinetics
Rates of enzyme-catalysed
reactions exhibit saturation:
Michaelis-Menten kinetics
Cooperativity: Hill function kinetics
Processes involving multiple
interacting components can
exhibit sigmoidal activity
(e.g. cooperative binding of O2
to hemoglobin)
Hill function (empirical fit)
Gene Regulatory Networks
Modelling constitutive gene expression
mRNA dynamics:
protein dynamics:
Modelling constitutive gene expression
If mRNA dynamics are fast (compared to protein dynamics):
Treat mRNA in ‘quasi-steady-state’
Reduced model only describes protein concentration :
Regulated gene expression
Constitutive expression
Repressed expression
Regulated gene expression
Constitutive expression
Activated expression
Modelling regulated expression
transcription factor
operator
Modelling regulated expression
transcription factor
operator
Regulation by multiple transcription
factors
transcription
factor B
transcription
factor A
operator
Distribution of states:
If A=B=P, with cooperative binding:
Natural Gene Regulatory
Networks
Autoregulating genes
Autorepressor: (regulation
enhances robustness and
response timing)
Autoactivator: (bistable ON/OFF
switching behaviour)
Gene switch: lac operon
Gene autoactivates in
response to lactose
Gene switch: lysis/lysogeny
decision in phage lambda
cI
cro
Double negative feedback
locks in one of two states
Oscillatory gene network: the
Goodwin oscillator
Delayed negative feedback
leads to sustained oscillations
Oscillatory gene network: circadian
rhythm generator
Delayed negative
feedback leads
to sustained
oscillations
Model: Goldbeter, 1996
Developmental gene networks
Segmentation in Drosophila
Endomesoderm specification in
purple sea urchin (Davidson, Bolouri
et al.)
Engineered Gene Circuits
The Collins Toggle Switch
Gardner, Cantor, and Collins, Nature, 2000
Double
repression
locks in one of
two possible
states.
Inducers allow
transitions
Collins toggle switch:
implementation
The Repressilator
Elowitz and Leibler, Nature, 2000
Three-step
repression ring
generates delayed
negative feedback:
sustained oscillations
progeny
single cell
Repressilator: Implementation
Improved oscillator design:
relaxation oscillator
Interplay of positive and negative feedback lead
to robust sustained rise-and-crash oscillations
Stricker/Hasty oscillator
Stricker et al., Nature, 2008
Interplay of
positive and
negative
feedback lead
to robust
sustained
oscillations
Stricker/Hasty oscillator
implementation
Lab I
Goals:
1) Simulate a differential equation model
in XPPAUT
2) Determine sensitivity coefficients for a
simple network model
3) Explore system behaviour in models of
gene regulatory networks: the Collins
toggle switch, the repressilator, the
Hasty/Stricker oscillator
Lab I
1) Simulate a differential equation model in XPPAUT
a) Open XPPAUT with file Lab1.ode or generate your own file, content is simply:
par k=1
x’=-k*x
init x=1
done
b) Select Initialconds|Go (I|G). Resize the window with Window/zoom|Fit (W|F). Then choose Initialconds|New
(I|N) and run a simulation with initial value of x set to 0.5.
c) Open the param window, change the value of k to 1. Re-run your simulations from x(0)=1 and x(0)=0.5. How has
the behaviour changed?
2) Determine sensitivity coefficients for a simple network model
a) Open XPPAUT with file Lab2.ode
b) Select Initialconds|Go (I|G). Use the Data window to view all four species concentrations. Select Graphic
stuff|Add curve (G|A) to add additional time-series to the plot.
c) Open the param window. Explore the effect of changing the parameter values. Consider the sensitivities of the
steady state of [A] with respect to (i) k1; (ii) k2; and (ii) k3.
3) Explore system behaviour in models of gene regulatory networks: the Collins toggle switch, the repressilator, the
Hasty/Stricker oscillator
a) Open XPPAUT with one of the files Lab_toggle.ode, Lab_repressilator.ode, or Lab_hastyosc.ode
b) Select Initialconds|Go (I|G). Verify the system’s desired behaviour: oscillations or bistability (for bistability, modify
the initial conditions).
c) Explore the effect of modifying the model parameters (param window) on the desired behaviour
Tools for analysis of
dynamic mathematical
models
The phase plane
The phase plane
Example network:
Time series:
species concentrations plotted
against time
Phase portrait:
species concentrations plotted
against one another
The phase plane
Time series: multiple simulations
(range of initial concentrations)
Phase portrait: multiple simulations
reveal system behaviour
The phase plane
Direction field
Nullclines
(turning points)
Intersection of nullclines: steady state.
Stability
Stability
Example: symmetric
antagonistic network:
Case I: Non-symmetric inhibition strengths:
unique long-time (steady-state) behaviour
Stability
Example: symmetric
antagonistic network:
Case II: Symmetric inhibition strengths:
two potential long-time (steady-state) behaviours
Stability
Nullclines intersect
three times.
Intermediate steadystate is unstable,
other two are stable.
Stability
A bistable system
Stability
monostability
bistability
saddle point
Linearized Stability Analysis
Model:
Evaluate Jacobian at steady state:
Determine eigenvalues of Jacobian:
If all eigenvalues have negative real part, then the
steady state is stable
Oscillations
Oscillatory behavior
Example: autocatalytic pathway:
Case I: weak autocatalysis:
damped oscillations (settling to steady state)
Oscillatory behavior
Example: autocatalytic pathway:
Case II: strong autocatalysis:
sustained (limit cycle) oscillations
Limit cycle
Bifurcation analysis
Bifurcation diagrams
The location of steady states depends on model parameters
Plot of steady-state concentration against parameter value:
continuation diagram
Bifurcation diagram: bistability
Parameter values at which the
number or stability of steady
states change are called
bifurcation points
nullclines vary with parameter values
Bistable
network
Bifurcation diagram: oscillations
Autocatalytic
network
Damped (transient) oscillation
persistent (limit cycle) oscillation
Model-based design
The Collins Toggle Switch
repression
Network:
repression
Goal:
bistability
Model:
Model-based
design
analysis
Conclusions:
need strong expression, highly nonlinear repression
The Repressilator
Network:
Goal:
sustained
oscillations
Model:
Model-based
design
analysis
Conclusions:
need low leak, strong non-linearity, short-lived proteins/long-lived mRNA
Alternative Modelling
Frameworks
Alternative modelling frameworks
Compartmental
modelling
Boolean
modelling
Stochastic
modelling
Spatial
Modelling
(PDEs)
Lab II
Goal:
Model-based design analysis of the Collins
toggle switch and the Hasty/Stricker
oscillator
Lab II
1)
a)
b)
c)
d)
e)
Model-based design analysis of the Collins toggle switch
Open XPPAUT with Lab_toggle.ode.
Generate a phase plane (Viewaxis|2D -- select P1 and P2 for the axes, set the window to [0,4]x[0,4]). Select
Initialconds|mIce (I|I) and generate multiple trajectories
Generate the nullclines: Nullclines|New (N|N)
Open the param window and change b to b=1. Reconstruct trajectories and nullclines. (Erase the previous phase
portrait.) How has the behaviour changed?
Explore a range of values for each parameter. What changes improve the robustness of the bistability? Which
changes eliminate it? What design choices can you recommend?
2) Model-based design analysis of the Hasty/Stricker oscillator relaxation
a) Open XPPAUT with Lab_hastyosc.ode
b) Repeat (b) and (c) above. (Set the window to X,Y [0,4]x[0,4]).
c) Open the param window and change alpha to alpha=5. Reconstruct trajectories and nullclines. (Erase the previous
phase portrait.) How has the behaviour changed?
d) Explore a range of values for each parameter. What changes eliminate the oscillations? What design choices can
you recommend?
3) Stability analysis
Determine the robustness of steady-state stability by determining the eigenvalues of the system Jacobian: select Sing.
pts.|Mouse, click near an equilibrium, print eigenvalues (negative real part signifies stability)
4) Bifurcation analysis
Generate a bifurcation diagram. Set parameters to their default values. Run a trajectory to steady state. Hit (I|L) a few
times to ensure steady state has been reached. Open AUTO: F|A. Set the parameter of interest to Par1 in the
Parameter window. Set the window size in the Axes window. Set the parameter value range (Par Min and Par Max) in
the Numerics window. Tap Run.
```