Transcript D - Grenoble

Qualitative Modeling and Simulation of Genetic
Regulatory Networks using Piecewise-Linear
Differential Equations
Hidde de Jong and Delphine Ropers
INRIA Rhône-Alpes
655 avenue de l’Europe, Montbonnot
38334 Saint Ismier Cedex
France
Email: {Hidde.de-Jong,Delphine.Ropers}@inrialpes.fr
Overview
1. Genetic regulatory networks
2. Models of genetic regulatory networks

nonlinear differential equations

linear differential equations

piecewise-linear differential equations
3. Qualitative modeling, simulation, and validation using
piecewise-linear differential equations
4. Genetic Network Analyzer (GNA)
2
Escherichia coli: model organism
 Enteric bacterium Escherichia coli has been most-studied
organism in biology
« All cell biologists have two cells of interest: the one they are studying and
Escherichia coli »
Schaechter and Neidhardt (1996), Escherichia coli and Salmonella, ASM Press, 4
107 bacteria
2 μm
4300 genes
3
Bacterial cell and proteins
 Proteins are building blocks of cell
Cell membrane, enzymes, gene expression, …
4
Variation in protein levels
 Protein levels in cell are adjusted to specific environmental
conditions
Peng, Shimizu (2003),
App. Microbiol. Biotechnol., 61:163-178
2D gels
DNA
microarrays
Western blots
Ali Azam et al. (1999), J. Bacteriol., 181(20):6361-6370
5
Synthesis and degradation of proteins
RNA polymerase
DNA
transcription
ribosome
effector molecule
modified protein
translation
mRNA
protein
post-translational
modification
degradation
protease
6
Regulation of synthesis and degradation
transcription
factor
RNA polymerase
DNA
RBS
modified protein
small RNA
mRNA
kinase
protease
ribosome
response regulator
7
Example: σS in E. coli
 σS (RpoS) is sigma factor in E. coli and other bacteria
Subunit of RNA polymerase which recognizes specific promoters
 σS is regulated on different levels:

Transcription: repression by CRP·cAMP

Translation: increase in efficiency by binding of small RNAs DsrA, RprA

Activity: increase in promoter affinity of RNAP with σS by binding of Crl

Degradation: RssB targets σS for degradation by ClpXP
Adapted from: Hengge-Aronis (2002),
Microbiol. Mol. Biol. Rev., 66(3):373-395
8
Genetic regulatory networks
 Control of protein synthesis and degradation gives rise to
genetic regulatory networks
Networks of genes, RNAs, proteins, metabolites, and their interactions
P
P
fis
P1 P2nlpD
gyrAB
GyrI
P1-P’1P2
cya
P
gyrI
σS
FIS
GyrAB
CYA
Supercoiling
Activation
TopA
P5 P1-P4
CRP
tRNA
rRNA
topA
Carbon starvation network in E. coli
rpoS
P1 P2
P1 P2
rrn
Stress
signal
RssB
crp
PA rssA PB rssB
9
Modeling of genetic regulatory networks
 Abundant knowledge on components and interactions of
genetic regulatory networks
 Currently no understanding of how global dynamics emerges
from local interactions between components
 Shift from structure to behavior of genetic regulatory networks
« functional genomics », « integrative biology », « systems biology », …
 Mathematical methods supported by computer tools allow
modeling and simulation of genetic regulatory networks:

precise and unambiguous description of network

understanding through computer experiments

new predictions
10
Model formalisms
 Many formalisms to model genetic regulatory networks
Graphs
Boolean equations
precision
abstraction
Differential equations
Stochastic master equations
de Jong (2002), J. Comput. Biol., 9(1): 69-105
 ODEs with implicit assumptions and additional simplifications:

Continuous and deterministic dynamics

Lumping together protein synthesis and degradation in single step
11
Cross-inhibition network
 Cross-inhibition network consists of two genes, each coding
for transcription regulator inhibiting expression of other gene
protein A
protein B
gene a
promoter a
promoter b
gene b
 Cross-inhibition network is example of positive feedback,
important for differentiation
Thomas and d’Ari (1990), Biological Feedback
12
Nonlinear model of cross-inhibition network
A
B
b
a
xa = a f (xb)  a xa
.
xb = b f (xa)  b xb
.
xa = concentration protein A
xb = concentration protein B
a, b > 0, production rate constants
a, b > 0, degradation rate constants
f (x )
1
f (x) =
0


n
 +x
n
n
, 
> 0 threshold
x
13
Phase-plane analysis
 Analysis of steady states in phase plane
.
xa
xa = 0
.
xa = 0 : xa =
.
.
xb = 0
0
xb = 0 : xb =
a
a f (xb)
b
f (xa)
b
xb
 Two stable and one unstable steady state. System will
converge to one of two stable steady states
 System displays hysteresis effect: transient perturbation may
cause irreversible switch to another steady state
14
Construction of cross inhibition network
 Construction of cross inhibition network in vivo
Gardner et al. (2000), Nature, 403(6786): 339-342
 Differential equation model of network
.
u=
α1
1+vβ
–u
.
v=
α2
1+u
–v
15
Experimental test of model
 Experimental test of mathematical model (bistability and
hysteresis)
Gardner et al. (2000), Nature, 403(6786): 339-342
16
Bifurcation analysis
 Analysis of bifurcations caused by changes in control
parameter
.
.
xa
xa
.
xa = 0
xa = 0
xa
.
.
xa = 0
.
xb = 0
0
xb = 0
xb = 0
xb
0
xb
0
xb
value of b
 Change in control parameter may cause an irreversible switch
to another steady state
17
Bacteriophage  infection of E. coli
 Response of E. coli to phage 
infection involves decision between
alternative developmental pathways:
lysis and lysogeny
Ptashne, A Genetic Switch, Cell Press,1992
18
Control of phage  fate decision
 Cross-inhibition feedback plays key role in establishment of
lysis or lysogeny, as well as in induction of lysis after DNA
damage
Santillán, Mackey (2004), Biophys. J., 86(1): 75-84
19
Simple model of phage  fate decision
 Differential equation model of cross-inhibition feedback network
involved in phage  fate decision
mRNA and protein, delays, thermodynamic description of gene regulation
Santillán, Mackey (2004), Biophys. J., 86(1): 75-84
20
Analysis of phage  model
 Bistability (lysis and lysogeny) only occurs for certain parameter
values
 Switch from lysis to lysogeny involves bifurcation from one
monostable regime to another, due to change in degradation
constant
Santillán, Mackey (2004), Biophys. J., 86(1): 75-84
21
Extended model of phage  infection
 Differential equation model of the extended network underlying
decision between lysis and lysogeny
McAdams, Shapiro (1995), Science, 269(5524): 650-656
22
Evaluation nonlinear differential equations
 Pro: reasonably accurate description of underlying molecular
interactions
 Contra: for more complex networks, difficult to analyze
mathematically, due to nonlinearities
 Pro: approximate solution can be obtained through numerical
simulation
 Contra: simulation techniques difficult to apply in practice, due
to lack of numerical values for parameters and initial conditions
23
Linear model of cross-inhibition network
A
B
b
a
xa = a f (xb)  a xa
.
xb = b f (xa)  b xb
.
xa = concentration protein A
xb = concentration protein B
a, b > 0, production rate constants
a, b > 0, degradation rate constants
f (x) = 1  x / (2 ) , 
1
> 0,
x  2
f (x )
0
2
x
24
Phase-plane analysis
 Analysis of steady states in phase plane
.
xa
xa = 0
.
xa = 0 : xa =
.
.
xb = 0
0
xb = 0 : xb =
a
a f (xb)
b
f (xa)
b
xb
 Single unstable steady state.
 Linear differential equations too simple to capture dynamic
phenomena of interest: no bistability and no hysteresis
25
Evaluation of linear differential equations
 Pro: analytical solution exists, thus facilitating analysis of
complex systems
 Contra: too simple to capture important dynamical
phenomena of regulatory network, due to neglect of nonlinear
character of interactions
26
Piecewise-linear model of cross-inhibition
A
B
a
xa = a f (xb)  a xa
.
xb = b f (xa)  b xb
.
b
xa = concentration protein A
xb = concentration protein B
a, b > 0, production rate constants
a, b > 0, degradation rate constants
f (x )
1
0
f (x) = s(x, ) =

1,
x<
0,
x>
x
Glass and Kauffman (1973), J. Theor. Biol., 39(1):103-129
27
PL models and gene regulatory logic
 Step function expressions correspond to Boolean functions
used to express gene regulatory logic
Thomas and d’Ari (1990), Biological Feedback
A
A
B
B
a
b
.
.
x
xa  a s-(xb , b ) – a xa
b
 b s-(xa , a) – b xb
condition gene a: (xb < b )
condition gene b: (xa < a )
.
.
x
b
a
xa  a s-(xa , a2) s-(xb , b ) – a xa
b
 b s-(xa , a1) – b xb
condition gene a: (xa < a2) 
condition gene b: (xa < a1)
(xb < b )
28
Phase-plane analysis
 Analysis of dynamics of PL models in phase space
κb/γb
κb/γb
xa
κa/γa
xa
a
a
M3
M1
b
0
M1:
.
.
x
xb
xa  a – a xa
b
 b – b xb
.
.
x
b
0
M3:
xb
.x  –  x
a
a a
.
xb  b – b xb
xa  a s-(xb , b ) – a xa
b
 b s-(xa , a) – b xb
29
Phase-plane analysis
 Analysis of dynamics of PL models in phase space
κb/γb
κb/γb
xa
κa/γa
xa
a
a
κa/γa
M5
M2
0
b
xb
.
.
x
0
b
xb
xa  a s-(xb , b ) – a xa
b
 b s-(xa , a) – b xb
 Extension of PL differential equations to differential inclusions
Gouzé, Sari (2002), Dyn. Syst., 17(4):299-316
using Filippov approach
30
Phase-plane analysis
 Global phase-plane analysis by combining analyses in local
regions of phase plane
.
xa
.
xa = 0
xa
xa = 0
a
.
.
xb = 0
xb = 0
0
b
xb
0
xb
 Piecewise-linear model good approximation of nonlinear
model, retaining properties of bistability and hysteresis
31
Qualitative analysis using PL models
 Hyperrectangular phase space partition: unique derivative
sign pattern in regions
 Qualitative abstraction yields state transition graph
Shift from continuous to discrete picture of network dynamics
D22
xa D19
D16
a D10
D23 D24
D25
D20 D21
D17 D18
D11 D12 D13D14 D15
D1
0
D2 D3 D4
D5
D6 D7 D8
D9
b
xb
D22
D23
D24
D19
D20
D21
D16
D10
D17
D11
.
D25
D1:
D18
D12
D13
D15
D14
xa > 0
.
xb > 0
.
x >0
D17: .a
xb < 0
.
D1
D2
D6
D5
D3 D4
D7
D9
D19:
xa = 0
.
xb = 0
D8
de Jong, Gouzé et al. (2004), Bull. Math. Biol., 66(2):301-340
32
Qualitative analysis using PL models
 Paths in state transition graph represent possible qualitative
behaviors
.
x >0
D1: .a
xb > 0
.
x >0
D17: .a
xb < 0
D22
D23
D24
D19
D20
D21
D17
D18
D16
D10
D11
κa/γa
D25
D12
D13
a
D15
D14
.
D19:
xa = 0
.
xb = 0
D1
D11 D17 D19
D1
D11 D17 D19
κb/γb
b
D1
D2
D6
D5
D3 D4
D7
D9
D8
33
Qualitative analysis using PL models
 State transition graph invariant for parameter constraints
κb/γb
xa
a
κa/γa
D11
D12
D11 D12
D1
D1
0
D3
b
D3
0 < a < a/a
0 < b < b/b
xb
34
Qualitative analysis using PL models
 State transition graph invariant for parameter constraints
κb/γb
xa
a
κa/γa
D11 D12
D11
D12
D1
D1
0
D3
b
D3
0 < a < a/a
0 < b < b/b
xb
35
Qualitative analysis using PL models
 State transition graph invariant for parameter constraints
κb/γb
xa
a
κa/γa
D11
D12
D11 D12
D1
D1
0
D3
0 < b < b/b
κb/γb
κa/γa
a
D3
0 < a < a/a
D11
D11
0 < a < a/a
D1
0 < b/b < b
D1
0
b
xb
36
Validation of qualitative models
 Predictions well adapted to comparison with available
experimental data: changes of derivative sign patterns
xa
0
time
xb
0
Concistency?
D22
D23
D24
D19
D20
D21
D17
D18
D16
D10
.
.
x >0
xa > 0
b
. time
.
x >0
xa < 0
b
D11
D25
D12
D1
D2
D6
D13
D15
D14
D5
D3 D4
D7
D9
D8
 Model validation: comparison of derivative sign patterns in
observed and predicted behaviors
 Need for automated and efficient tools for model validation
37
Validation of qualitative models
 Predictions well adapted to comparison with available
experimental data: changes of derivative sign patterns
xa
0
time
xb
0
.
x >0
.
x >0
a
b
. time
x <0
.
x >0
a
b
D22
D23
D24
D19
D20
D21
D17
D18
Concistency?
D16
Yes
D10
D11
.
D25
D12
D13
x >0
D1: .a
xb > 0
D15
D14
.
D17:
xa > 0
.
xb < 0
D19:
xa = 0
.
xb = 0
.
D1
D2
D6
D5
D3 D4
D7
D9
D8
 Model validation: comparison of derivative sign patterns in
observed and predicted behaviors
 Need for automated and efficient tools for model validation
38
Model-checking approach
 Dynamic properties of system can be expressed in temporal
logic (CTL)
x
a
.
.
There Exists a Future state where xa > 0 and xb > 0
and starting from that state,
.
.
there Exists a Future state where xa < 0 and xb > 0
.
.
.
.
EF(xa > 0  xb > 0  EF(xa < 0  xb > 0) )
0
time
xb
0
.
.
x >0
xa > 0
b
. time
.
x >0
xa < 0
b
 Model checking is automated technique for verifying that state
transition graph satisfies temporal-logic statements
 Computer tools are available to perform efficient and reliable
model checking (NuSMV, CADP, …)
39
Validation using model checking
 Compute state transition graph using qualitative simulation
xa
D23
D24
D19
D20
D21
D17
D18
D16
0
time
D10
xb
0
D22
.
.
x >0
xa > 0
b
. time
.
x >0
xa < 0
Concistency?
D11
D25
D12
D1
D2
D6
D13
D15
D14
D5
D3 D4
D7
D9
D8
b
 Use of model checkers to verify whether experimental data and
predictions are consistent
Batt et al. (2005), Bioinformatics, 21(supp. 1): i19-i28
40
Validation using model checking
 Compute state transition graph using qualitative simulation
.
.
.
.
EF(xa > 0  xb > 0  EF(xa < 0  xb > 0) )
D22
D23
D24
D19
D20
D21
D17
D18
D16
D10
Concistency?
Yes
Model corroborated
D11
D25
D12
D1
D2
D6
D13
D15
D14
D5
D3 D4
D7
D9
D8
 Use of model checkers to verify whether experimental data and
predictions are consistent
Batt et al. (2005), Bioinformatics, 21(supp. 1): i19-i28
41
Analysis of attractors of PL systems
 Search of steady states of PL systems in phase space
xa
D22
D23
D24
D19
D20
D21
D17
D18
D16
a
D10
D11
D12
D1
0
b
xb
D25
D2
D6
D13
D15
D14
D3 D4
D7
D5
D9
D8
42
Analysis of attractors of PL systems
 Search of steady states of PL systems in phase space
xa
D22
D23
D24
D19
D20
D21
D17
D18
D16
a
D10
D11
D12
D1
0
b
xb
D25
D2
D6
D13
D15
D14
D3 D4
D7
D5
D9
D8
 Analysis of stability of steady states, using local properties of
state transition graph
Casey et al. (2006), J. Math Biol., 52(1):27-56
Definition of stability of equilibrium points on surfaces of discontinuity
43
Genetic Network Analyzer (GNA)
 Qualitative simulation method implemented in Java: Genetic
Network Analyzer (GNA)
Distribution by
Genostar SA
de Jong et al. (2003),
Bioinformatics, 19(3):336-344
Batt et al. (2005), Bioinformatics,
21(supp. 1): i19-i28
http://www-helix.inrialpes.fr/gna
44
Applications of GNA
 Qualitative simulation method used to analyze various
bacterial regulatory networks:

initiation of sporulation in Bacillus subtilis
de Jong, Geiselmann et al. (2004), Bull. Math. Biol., 66(2):261-300

quorum sensing in Pseudomonas aeruginosa
Viretta and Fussenegger, Biotechnol. Prog., 2004, 20(3):670-678

carbon starvation response in Escherichia coli
Ropers et al., Biosystems, 2006, 84(2):124-152

onset of virulence in Erwinia chrysanthemi
Sepulchre et al., J. Theor. Biol., 2006, in press
45
Evaluation of PL differential equations
 Pro: captures important dynamical phenomena of network, by
suitable approximation of nonlinearities
 Pro: qualitative analysis of dynamics possible, due to favorable
mathematical properties
 Contra: restricted class of models, not directly applicable to
type of functions found in, for example, metabolism
46
Contributors and sponsors
Grégory Batt, Boston University, USA
Hidde de Jong, INRIA Rhône-Alpes, France
Hans Geiselmann, Université Joseph Fourier, Grenoble, France
Jean-Luc Gouzé, INRIA Sophia-Antipolis, France
Radu Mateescu, INRIA Rhône-Alpes, France
Michel Page, INRIA Rhône-Alpes/Université Pierre Mendès France, Grenoble, France
Corinne Pinel, Université Joseph Fourier, Grenoble, France
Delphine Ropers, INRIA Rhône-Alpes, France
Tewfik Sari, Université de Haute Alsace, Mulhouse, France
Dominique Schneider, Université Joseph Fourier, Grenoble, France
Ministère de la Recherche,
IMPBIO program
European Commission,
FP6, NEST program
INRIA, ARC program
Agence Nationale de la
Recherche, BioSys program
47