BioNetGen Tutorial

Download Report

Transcript BioNetGen Tutorial

A BIONETGEN Tutorial
Outline
• Downloading BNG
• Installing BNG
• Running BNG
– text-based interface
– graphical interface – RULEBUILDER
– web-based interface – GETBONNIE
• Simple example
• Extending the example
• Detailed model of TLR4 signaling
Creating a BNG model
6. Create a BioNetGen model (with .bngl extension) in a text file using
your favorite editor.
Running BNG
7. Open a terminal window, navigate to the file location, and type path
to BNG2.pl followed by filename.
The BNG logfile
8. (optional) Piping output to a logfile.
~/BioNetGen2/Perl2/BNG2.pl example1.bngl > & example1.log
BioNetGen version 2.0.48
readFile::Reading from file example1.bngl
Read 11 parameters.
Read 3 species.
Adding P as allowed state of component Y of molecule R
Read 3 observable(s).
Read 5 reaction rule(s).
Iteration 0: 3 species
0 rxns 0.00e+00 CPU s
Iteration 1: 4 species
1 rxns 1.00e-02 CPU s
Iteration 2: 5 species
3 rxns 0.00e+00 CPU s
Iteration 3: 6 species
5 rxns 1.00e-02 CPU s
Iteration 4: 9 species
9 rxns 1.00e-02 CPU s
Iteration 5: 12 species 20 rxns 3.00e-02 CPU s
Iteration 6: 14 species 32 rxns 3.00e-02 CPU s
Iteration 7: 14 species 36 rxns 2.00e-02 CPU s
Cumulative CPU time for each rule
Rule 1: 6 reactions 3.00e-02 CPU s 5.00e-03 CPU s/rxn
Rule 2: 12 reactions 5.00e-02 CPU s 4.17e-03 CPU s/rxn
Rule 3: 3 reactions 0.00e+00 CPU s 0.00e+00 CPU s/rxn
Rule 4: 5 reactions 2.00e-02 CPU s 4.00e-03 CPU s/rxn
Rule 5: 10 reactions 1.00e-02 CPU s 1.00e-03 CPU s/rxn
Total : 36 reactions 1.10e-01 CPU s 3.06e-03 CPU s/rxn
Wrote network to example1.net.
CPU TIME: generate_network 0.1 s.
Network simulation using ODEs
Running run_network on phoebe
Overview of inputs and outputs
NETWORK
.cdat
ODE
BNGL
BIONETGEN
.net
SSA
.gdat
.log
Elements of a BNGL file
•
•
•
•
•
•
parameters
molecule types
seed species
reaction rules
observables
actions
A simple example:
Ligand-receptor binding
EGF
Vo
V
EGFR
A simple example:
parameters
begin parameters
# Physical and geometric constants
NA 6.0e23 #Avogadro’s num
f
1 #scaling factor
Vo
f*1e-10 # L
V
f*3e-12 # L
EGF
Vo
V
EGFR
# Initial concentrations
EGF0 2e-9*NA*Vo #nM
EGFR0 f*1.8e5 #copy per cell
# Rate constants
kp1 9.0e7/(NA*Vo) #input /M/sec
km1 0.06
#/sec
end parameters
A simple example:
parameters
Summary
EGF
Vo
V
EGFR
Concentration units : copies per cell
• Multiply concentrations by NA*Vrxn
• Don’t scale first order rate constants
• Divide second order rate constants by NA*Vrxn
Following this recipe allows switching between
ODE and stochastic simulation without
parameter modification.
A simple example:
molecule types
R
EGF
begin molecule types
EGF(R)
EGFR(L,CR1,Y1068~U~P)
L
CR1
Y1068
EGFR
U
P
end molecule types
A simple example:
seed species
R
EGF
begin seed species
L
CR1
Y1068~U
EGFR
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
Pattern basics
Basic Definition: A set of molecules that may be partially
specified and selects a set of species through its allowed
mappings.
Pattern
EGF(R)
EGFR(L)
EGF(R)
EGFR(L,CR1,Y1068~U)
Mapping
Species
EGFR(L,CR1,Y1068~P)
A simple example:
reaction rules
+
EGF
EGFR
k+1
k-1
EGF
EGFR
begin reaction rules
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
end reaction rules
A simple example:
actions
begin parameters
…
end parameters
begin molecule types
EGF(R)
EGFR(L,CR1,Y1068~U~P)
end molecule types
begin seed species
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
begin reaction rules
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
end reaction rules
# actions
generate_network({overwrite=>1});
Generates network of species and reactions by iterative application of rules
starting with seed species
A simple example:
Application of the forward binding rule
Reaction rule
EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
Add bond
Step 1: Match patterns onto reactant species.
Step 2: Copy selected species to product side.
Step 3: Apply transformation specified by rule.
A simple example:
Application of the forward binding rule
Reaction rule
EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
Add bond
Step 1: Match patterns onto reactant species.
EGF(R) + EGFR(L,CR1,Y1068~U)
A simple example:
Application of the forward binding rule
Reaction rule
EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
Add bond
Step 1: Match patterns onto reactant species.
EGF(R) + EGFR(L,CR1,Y1068~U)
Step 2: Copy selected species to product side.
-> EGF(R) + EGFR(L,CR1,Y1068~U)
A simple example:
Application of the forward binding rule
Reaction rule
EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
Add bond
Step 1: Match patterns onto reactant species.
EGF(R) + EGFR(L,CR1,Y1068~U)
Step 2: Copy selected species to product side.
-> EGF(R) + EGFR(L,CR1,Y1068~U)
Step 3: Apply transformation specified by rule.
-> EGF(R!1).EGFR(L!1,CR1,Y1068~U) kp1
New species
A simple example:
Application of the forward binding rule
Reaction rule
EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
Add bond
EGF(R) + EGFR(L,CR1,Y1068~U)
-> EGF(R!1).EGFR(L!1,CR1,Y1068~U) kp1
New species
New reaction
A simple example:
Application of the reverse binding rule
Reaction rule
EGF(R!1).EGFR(L!1) -> EGF(R) + EGFR(L) km1
Delete bond
EGF(R!1).EGFR(L!1,CR1,Y1068~U)
A simple example:
Application of the reverse binding rule
Reaction rule
EGF(R!1).EGFR(L!1) -> EGF(R) + EGFR(L) km1
Delete bond
EGF(R!1).EGFR(L!1,CR1,Y1068~U)
-> EGF(R!1).EGFR(L!1,CR1,Y1068~U)
A simple example:
Application of the reverse binding rule
Reaction rule
EGF(R!1).EGFR(L!1) -> EGF(R) + EGFR(L) km1
Delete bond
EGF(R!1).EGFR(L!1,CR1,Y1068~U)
-> EGF(R) + EGFR(L,CR1,Y1068~U) km1
New reaction2
A simple example:
Overview of generate_network
seed species
EGF(R)
EGFR(L,CR1,Y1068~U)
Iteration 1:
rule1f: EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
A simple example:
Overview of generate_network
seed species
EGF(R)
EGFR(L,CR1,Y1068~U)
Iteration 1:
rule1f: EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
new reaction: EGF(R) + EGFR(L,CR1,Y1068~U) ->
EGF(R!1).EGFR(L!1,CR1,Y1068~U) kp1
new species: EGF(R!1).EGFR(L!1,CR1,Y1068~U)
A simple example:
Overview of generate_network
seed species
EGF(R)
EGFR(L,CR1,Y1068~U)
Iteration 1:
rule1f: EGF(R) + EGFR(L) -> EGF(R!1).EGFR(L!1) kp1
new reaction: EGF(R) + EGFR(L,CR1,Y1068~U) ->
EGF(R!1).EGFR(L!1,CR1,Y1068~U) kp1
new species: EGF(R!1).EGFR(L!1,CR1,Y1068~U)
Network after Iteration 1
species
1 EGF(R)
2 EGFR(L,CR1,Y1068~U)
3 EGF(R!1).EGFR(L!1,CR1,Y1068~U)
reactions
1 1,2 3 kp1
A simple example:
Overview of generate_network
From iteration 1
species
1 EGF(R)
2 EGFR(L,CR1,Y1068~U)
3 EGF(R!1).EGFR(L!1,CR1,Y1068~U)
Iteration 2:
rule1r: EGF(R!1).EGFR(L!1) -> EGF(R) + EGFR(L) km1
new reaction: EGF(R!1).EGFR(L!1,CR1,Y1068~U)
-> EGF(R) + EGFR(L,CR1,Y1068~U) km1
Network after Iteration 2
species
1 EGF(R)
2 EGFR(L,CR1,Y1068~U)
3 EGF(R!1).EGFR(L!1,CR1,Y1068~U)
reactions
1 1,2 3 kp1 # rule1f
2 3 1,2 km1 # rule1r
A simple example:
Overview of generate_network
From iteration 1
species
1 EGF(R)
2 EGFR(L,CR1,Y1068~U)
3 EGF(R!1).EGFR(L!1,CR1,Y1068~U)
Iteration 2:
rule1r: EGF(R!1).EGFR(L!1) -> EGF(R) + EGFR(L) km1
new reaction: EGF(R!1).EGFR(L!1,CR1,Y1068~U)
-> EGF(R) + EGFR(L,CR1,Y1068~U) km1
Network after Iteration 2 - final
species
1 EGF(R)
2 EGFR(L,CR1,Y1068~U)
3 EGF(R!1).EGFR(L!1,CR1,Y1068~U)
reactions
1 1,2 3 kp1 # rule1f
2 3 1,2 km1 # rule1r
A simple example:
Running BNG2.pl
BNGL
BNG2.pl
.log
generate_network
.net
A simple example:
Running BNG2.pl
phoebe% ~/BioNetGen2/Perl2/BNG2.pl SimpleExample.bngl
/Users/faeder/BioNetGen2/Perl2/BNG2.pl
BioNetGen version 2.0.48
readFile::Reading from file SimpleExample.bngl
Read 8 parameters.
Read 2 molecule types.
Read 2 species.
Read 1 reaction rule(s).
WARNING: Removing old network file SimpleExample.net.
Iteration 0: 2 species
0 rxns 0.00e+00 CPU s
Iteration 1: 3 species
1 rxns 0.00e+00 CPU s
Iteration 2: 3 species
2 rxns 0.00e+00 CPU s
Cumulative CPU time for each rule
Rule 1: 2 reactions 0.00e+00 CPU s 0.00e+00 CPU s/rxn
Total : 2 reactions 0.00e+00 CPU s 0.00e+00 CPU s/rxn
Wrote network to SimpleExample.net.
A simple example:
The .net file
The .net
file
begin molecule types
1 EGF(R)
2 EGFR(L,CR1,Y1068~U~P)
end molecule types
begin species
1 EGF(R)
EGF0
2 EGFR(CR1,L,Y1068~U)
EGFR0
3 EGF(R!1).EGFR(CR1,L!1,Y1068~U) 0
end species
begin reaction rules
Rule1: \
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
# Bind(0.0.0,0.1.0)
# Reverse
# Unbind(0.0.0,0.1.0)
end reaction rules
begin reactions
1 1,2 3 kp1 #Rule1
2 3 1,2 km1 #Rule1r
end reactions
SimpleExample.net
Translation of the network model into
ordinary differential equations
Dependent variables
species
1 EGF(R)
EGF0
2 EGFR(CR1,L,Y1068~U)
EGFR0
3 EGF(R!1).EGFR(CR1,L!1,Y1068~U) 0
reactions
1 1,2 3 kp1 #Rule1
2 3 1,2 km1 #Rule1r
Flux terms
f1
f2
Exporting the model
BNGL
BNG2.pl
.log
generate_network
.net
writeformat
M-file
SBML
Exporting model equations to Matlab:
writeMfile
begin parameters
…
end parameters
begin molecule types
EGF(R)
EGFR(L,CR1,Y1068~U~P)
end molecule types
begin seed species
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
begin reaction rules
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
end reaction rules
# actions
generate_network({overwrite=>1});
writeMfile();
Ordinary differential equations in the
M-file
SimpleExample.m
% M-file for model SimpleExample created by
BioNetGen 2.0.48
function [t_out,obs_out,x_out]=
SimpleExample(tend)
Nspecies=3;
Nreactions=2;
obs_out=zeros(1);
% Parameters
NA=6.0e23;
f=1;
Vo=f*1e-10;
V=f*3e-12;
EGF0=(2e-9*NA)*Vo;
EGFR0=f*1.8e5;
kp1=3.0e6/(NA*V);
km1=0.06;
% Intial concentrations
x0= [ EGF0; EGFR0; 0;];
% Reaction flux function
function f= flux(t,x)
f(1)=kp1*x(1)*x(2);
f(2)=km1*x(3);
end
% Derivative function
function d= xdot(t,x)
f=flux(t,x);
d(1,1)= -f(1) +f(2);
d(2,1)= -f(1) +f(2);
d(3,1)= +f(1) -f(2);
end
snames={'EGF(R)','EGFR(CR1,L,Y1068~U)','EGF(R!1).
EGFR(CR1,L!1,Y1068~U)'};
% Integrate ODEs
[t_out,x_out]= ode15s(@xdot, [0 tend], x0);
% plot species concentrations
plot(t_out,x_out);
legend(snames);
end
Running exported model as a Matlab
function
>> SimpleExample(120)
Export to the Systems Biology Markup
Language (SBML)
• SBML is an XML specification for biological
models
• Many simulation tools read and write SBML
• Supports standard reaction network models
• generate_network must be invoked prior to
export
• Support for rule-based models will appear in
Level 3 (current level is 2+)
• Website: http://sbml.org
Running the simulation within
BioNetGen
NETWORK
generate_network
BNGL
BNG2.pl
simulate_ode
.cdat
ODE
.net
SSA
simulate_ssa
.log
.gdat
Running the simulation within
BioNetGen
begin parameters
…
end parameters
begin molecule types
EGF(R)
EGFR(L,CR1,Y1068~U~P)
end molecule types
begin seed species
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
begin reaction rules
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
end reaction rules
# actions
generate_network({overwrite=>1});
saveConcentrations();
simulate_ode({suffix=>ode,t_end=>120,n_steps=>50});
Running the simulation within
BioNetGen
begin parameters
…
end parameters
begin molecule types
EGF(R)
EGFR(L,CR1,Y1068~U~P)
end molecule types
begin seed species
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
BNG, like many other systems biology
begin reaction rules
applications, uses the CVODE library for
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
solving ODE’s.
end reaction rules
# actions
generate_network({overwrite=>1});
saveConcentrations();
Sparse option can be used to solve
systems up to about 50,000 equations.
simulate_ode({suffix=>ode,t_end=>120,n_steps=>50});
The .cdat file
SimpleExample_ode.cdat
#
time
0.0000000000000000e+00
2.3999999999999999e+00
4.7999999999999998e+00
7.1999999999999993e+00
1
2
1.2000000000000000e+05
6.9387285118773594e+04
5.0363939566816996e+04
4.1819631802026437e+04
3
1.8000000000000000e+05
1.2938728511877365e+05
1.1036393956681706e+05
1.0181963180202652e+05
0.0000000000000000e+00
5.0612714881226493e+04
6.9636060433183069e+04
7.8180368197973614e+04
Your favorite plotting package
Some useful options
• PhiBPlot – Java-based (JFreeChart) plotting utility provided with BNG
• xmgrace – Unix / X windows graphing package great for XY plots
• Matlab, Mathematica, gnuplot, Excel, etc.
The .cdat file in PhiBPlot
The .cdat file in PhiBPlot
These
numbers
correspond to
indices of
species in .net
file
The .cdat file in PhiBPlot
These
numbers
correspond to
indices of
species in .net
file
begin species
1 EGF(R)
EGF0
2 EGFR(CR1,L,Y1068~U)
EGFR0
3 EGF(R!1).EGFR(CR1,L!1,Y1068~U) 0
end species
Stochastic simulations in BioNetGen
• Network simulation engine uses Gillespie
Direct method
• Adequate performance for networks with up
to ~50,000 species.
• May be used in ‘on-the-fly’ mode, which
generates species and reactions as needed
Stochastic simulations in BioNetGen
• Network simulation engine uses Gillespie
Direct method
• Adequate performance for networks with up
to ~50,000 species.
• May be used in ‘on-the-fly’ mode, which
generates species and reactions as needed
Comparison between ODE and
stochastic simulations
begin parameters
…
end parameters
begin molecule types
EGF(R)
EGFR(L,CR1,Y1068~U~P)
end molecule types
begin seed species
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
begin reaction rules
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
end reaction rules
# actions
generate_network({overwrite=>1});
saveConcentrations();
simulate_ode({suffix=>ode,t_end=>120,n_steps=>50});
resetConcentrations();
simulate_ssa({suffix=>ssa,t_end=>120,n_steps=>50});
Viewing comparison in PhiBPlot
Viewing comparison in PhiBPlot
Viewing comparison in PhiBPlot
Noise is minimal
because species
concentrations are
high.
Viewing comparison in PhiBPlot
begin parameters
# Physical and geometric constants
NA 6.0e23 #Avogadro’s num
f
1 #scaling factor
Vo
f*1e-10 # L
V
f*3e-12 # L
Viewing comparison in PhiBPlot
begin parameters
# Physical and geometric constants
NA 6.0e23 #Avogadro’s num
f
0.001 #scaling factor
Vo
f*1e-10 # L
V
f*3e-12 # L
Viewing comparison in PhiBPlot
Let’s add another rule: dimerization
+
EGF
EGFR
k+2
k-2
Let’s add another rule: dimerization
+
EGF
k+2
k-2
EGFR
EGF(R!2).EGFR(L!2,CR1) + EGFR(L!3,CR1).EGF(R!3)
Selects EGFR with free dimerization domain (CR1) and
bound ligand-binding domain (L).
Let’s add another rule: dimerization
+
EGF
k+2
k-2
EGFR
EGF(R!2).EGFR(L!2,CR1) + EGFR(L!3,CR1).EGF(R!3)
A shorthand way to select the same properties uses the
bond wildcard (‘!+’) and omits the EGF binding partner.
EGFR(L!+,CR1) + EGFR(L!+,CR1)
These patterns select EGFR molecules with free CR1 sites
and their L sites bound to anything.
Let’s add another rule: dimerization
+
EGF
k+2
k-2
EGFR
Full rule
EGFR(L!+,CR1) + EGFR(L!+,CR1) <-> EGFR(L!+,CR1!1).EGFR(L!+,CR1!1) \
kp2,km2
Let’s add another rule: dimerization
+
EGF
k+2
k-2
EGFR
Full rule
EGFR(L!+,CR1) + EGFR(L!+,CR1) <-> EGFR(L!+,CR1!1).EGFR(L!+,CR1!1) \
kp2,km2
Need to add kp2 and km2 to parameters block.
Application of the dimerization rule
EGFR(L!+,CR1) + EGFR(L!+,CR1)
EGFR(L!1,CR1,Y1068~U).
EGF(R!1)
EGFR(L!1,CR1,Y1068~U).
EGF(R!1)
add bond
EGF(R!1).EGFR(L!1,CR1!3,Y1068~U).EGFR(L!2,CR1!3,Y1068~U).EGF(R!2)
species 4
New reaction
reactions
…
3 3,3 4 0.5*kp2
Application of the dimerization rule
EGFR(L!+,CR1) + EGFR(L!+,CR1)
EGFR(L!1,CR1,Y1068~U).
EGF(R!1)
EGFR(L!1,CR1,Y1068~U).
EGF(R!1)
add bond
EGF(R!1).EGFR(L!1,CR1!3,Y1068~U).EGFR(L!2,CR1!3,Y1068~U).EGF(R!2)
species 4
New reaction
reactions
…
3 3,3 4 0.5*kp2
Symmetry factor
An unexpected result
.bngl file
…
begin seed species
EGF(R)
EGF0
EGFR(L,CR1,Y1068~U) EGFR0
end seed species
begin reaction rules
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
EGFR(L!+,CR1) + EGFR(L!+,CR1) <-> \
EGFR(L!+,CR1!1).EGFR(L!+,CR1!1) kp2, km2
end reaction rules
generate_network({overwrite=>1});
logfile
…
Iteration
Iteration
Iteration
Iteration
Iteration
Iteration
…
0:
1:
2:
3:
4:
5:
2 species
3 species
4 species
5 species
6 species
6 species
0 rxns
1 rxns
3 rxns
5 rxns
7 rxns
8 rxns
0.00e+00 CPU s
0.00e+00 CPU s
0.00e+00 CPU s
1.00e-02 CPU s
1.00e-02 CPU s
0.00e+00 CPU s
Expect network
generation to stop here!
(+ one more rxn)
Checking species in the .net file
begin species
1 EGF(R)
EGF0
2 EGFR(CR1,L,Y1068~U)
EGFR0
3 EGF(R!1).EGFR(CR1,L!1,Y1068~U)
0
4 EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U) 0
5 EGF(R!1).EGFR(CR1!2,L!1,Y1068~U).EGFR(CR1!2,L,Y1068~U)
0
6 EGFR(CR1!1,L,Y1068~U).EGFR(CR1!1,L,Y1068~U)
0
end species
Dimers with 0 or 1 EGF molecule arise because EGF is allowed to
dissociate whether or not EGFR is dimerized.
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
state of CR1 domain is not specified
Checking species in the .net file
begin species
1 EGF(R)
EGF0
2 EGFR(CR1,L,Y1068~U)
EGFR0
3 EGF(R!1).EGFR(CR1,L!1,Y1068~U)
0
4 EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U) 0
5 EGF(R!1).EGFR(CR1!2,L!1,Y1068~U).EGFR(CR1!2,L,Y1068~U)
0
6 EGFR(CR1!1,L,Y1068~U).EGFR(CR1!1,L,Y1068~U)
0
end species
Dimers with 0 or 1 EGF molecule arise because EGF is allowed to
dissociate whether or not EGFR is dimerized.
EGF(R) + EGFR(L) <-> EGF(R!1).EGFR(L!1) kp1, km1
state of CR1 domain is not specified
“Expected” behavior can be recovered by modifying the rule:
EGF(R) + EGFR(L,CR1) <-> EGF(R!1).EGFR(L!1,CR1) kp1, km1
Only EGFR with free CR1 will bind and release ligand.
Creating an observable:
Receptors in dimers
Like rules, Observables select species based on pattern
matching.
EGFR(CR1!+)
matches receptors in dimers
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U)
The concentrations of the matching species are stored in a
variable with a specified name, computed by summing the
concentrations of matching species.
begin observables
Molecules Rdim EGFR(CR1!+)
end observables
Creating an observable:
Receptors in dimers
Like rules, Observables select species based on pattern
matching.
EGFR(CR1!+)
matches receptors in dimers
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U)
begin observables
Molecules Rdim EGFR(CR1!+)
end observables
generate_network
.net file
begin groups
1 Rdim
end groups
number of matches
2*4
species index
The .gdat file
#
time
0.00000000e+00
2.40000000e+00
4.80000000e+00
…
Rdim
0.00000000e+00
1.32882737e+00
5.16753986e+00
The .gdat file
Rdim
2x
Dimers
(Species 4)
A state change:
Transphosphorylation
Kinase domain in
cytoplasmic tail of
dimerized receptor
phosphorylates tyrosines
of trans receptor.
P
EGFR(CR1!+,Y1068~U) -> EGFR(CR1!+,Y1068~P) kp3
State change
A state change:
Transphosphorylation
Kinase domain in
cytoplasmic tail of
dimerized receptor
phosphorylates tyrosines
of trans receptor.
P
EGFR(CR1!+,Y1068~U) -> EGFR(CR1!+,Y1068~P) kp3
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U)
State change
Because of symmetry, produce two copies of same reaction.
Reaction multiplicity
Transphosphorylation rule generates
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U)
Species 4
kp3
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~P).EGFR(CR1!3,L!2,Y1068~U)
Species 5
and
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~U)
Species 4
kp3
EGF(R!1).EGF(R!2).EGFR(CR1!3,L!1,Y1068~U).EGFR(CR1!3,L!2,Y1068~P)
Species 5
which are combined into
reactions
…
5 4 5 2*kp3 #Rule3
Dephosphorylation
Phosphatase excess gives
rise to first order
dephosphorylation
process.
P
U
EGFR(Y1068~P) -> EGFR(Y1068~U) km3
Note the convention that a component
without a specified bonding state is
required to be unbound.
Dephosphorylation can occur only if
the site is unbound.
Binding of another protein to the site
protects it from dephosphorylation.
Monitoring site phosphorylation
Correct way to specify observable for receptor phosphorylation
observables
Molecules RP EGFR(Y1068~P!?)
selects
site in
state ‘P’
wildcard
allows site
to be
bound or
unbound
Monitoring site phosphorylation
Correct way to specify observable for receptor phosphorylation
observables
Molecules RP EGFR(Y1068~P!?)
Incorrect way to specify observable for receptor phosphorylation
observables
Molecules RP EGFR(Y1068~P)
Only monitors sites
that are
phosphorylated and
unbound
Adapter and effector binding
P
EGFR
5
Grb2
6
Sos1
cytosolic proteins
recruited to the
membrane by
receptor
phosphorylation
Additional rules
5 EGFR(Y1068~P) + Grb2(SH2) <-> EGFR(Y1068~P!1).Grb2(SH2!1) kp4,km4
6 Grb2(SH3) + Sos1(PxxP) <-> Grb2(SH3!1).Sos1(PxxP!1) kp5,km5
Observable
Molecules Sos1_act Sos1(PxxP!+)
Selects Sos1 bound to the
membrane
Adapter and effector binding
P
EGFR
5
Grb2
6
Sos1
cytosolic proteins
recruited to the
membrane by
receptor
phosphorylation
Additional rules
5 EGFR(Y1068~P) + Grb2(SH2) <-> EGFR(Y1068~P!1).Grb2(SH2!1) kp4,km4
6 Grb2(SH3) + Sos1(PxxP) <-> Grb2(SH3!1).Sos1(PxxP!1) kp5,km5
Observable
Molecules Sos1_act Sos1(PxxP!+)
Selects Sos1 bound to the
membrane - WRONG!!
Adapter and effector binding
P
EGFR
5
Grb2
6
Sos1
cytosolic proteins
recruited to the
membrane by
receptor
phosphorylation
Additional rules
5 EGFR(Y1068~P) + Grb2(SH2) <-> EGFR(Y1068~P!1).Grb2(SH2!1) kp4,km4
6 Grb2(SH3) + Sos1(PxxP) <-> Grb2(SH3!1).Sos1(PxxP!1) kp5,km5
Observable
Molecules Sos1_act Sos1(PxxP!1).Grb2(SH3!1,SH2!2).EGFR(Y1068!2)
Correct observable traces connectivity to membrane
bound molecule
A degradation reaction
begin molecule types
…
Null()
…
begin seed species
…
$Null
Simple receptor degradation
begin reaction rules
…
EGFR() -> Null kdeg
Deletes any species containing EGFR molecule,
including any attached molecules.
A degradation reaction
begin molecule types
…
Null()
…
begin seed species
…
$Null
Selective degradation of ligand and receptor
begin reaction rules
…
EGF(R!1).EGF(R!2).EGFR(L!1,CR1!3).EGFR(L!2,CR1!3) -> Trash() kdeg \
DeleteMolecules
This rule causes just the ligand-bound receptor dimer to be
degraded, recycling any components that are bound to the
cytoplasmic domains of the receptor.
Running the full model
EGF
1
1’
2
4 U
P
5
EGFR
3
Grb2
generate_network
…
Iteration
Iteration
Iteration
Iteration
Iteration
Iteration
Iteration
Iteration
0:
1:
2:
3:
4:
5:
6:
7:
5 species
7 species
8 species
9 species
13 species
18 species
23 species
23 species
0 rxns 0.00e+00 CPU s
2 rxns 0.00e+00 CPU s
5 rxns 1.00e-02 CPU s
8 rxns 1.00e-02 CPU s
14 rxns 2.00e-02 CPU s
35 rxns 6.00e-02 CPU s
66 rxns 9.00e-02 CPU s
86 rxns 9.00e-02 CPU s …
6
Sos1
Results
Equilibrating Grb2-Sos1 binding prior
to ligand addition
Grb2-Sos1 interaction
is constitutive
Grb2
6
Sos1
Equilibrating Grb2-Sos1 binding prior
to ligand addition
Grb2-Sos1 interaction
is constitutive
Grb2
6
Sos1
Combination of actions can be used to pre-equilibrate prior
to ligand addition
setConcentration("EGF(R)",0);
simulate_ode({suffix=>equil,t_end=>10000,n_steps=>50, \
steady_state=>1,sparse=>1});
setConcentration("EGF(R)",EGF0);
Complete summary of actions
Action/parameter
Type
generate_network
Description
Default
Generate species and reactions through iterative application of
rules to seed species
max_agg
int
Maximum number of molecules in one species
1e99
max_iter
int
Maximum number of iterations of rule
application
100
max_stoich
hash
Maximum number of molecules of specified
type in one species
-
overwrite
0/1
Overwrite existing .net file
0 (off)
print_iter
0/1
Print .net file after each iteration
0
prefixa
string
Set basenamec of .net file to string
basename of .bngl
file
suffixa
string
Append _string to basename of .net file
-
simulate_ode / simulate_ssa
Simulate current model/network
t_end
float
End time for simulation
required
t_start
float
Start time for simulation
0
n_steps
int
Number of times after t=0 at which to report
concentrations/observables
1
sample_times
array
Times at which to report
concentrations/observables (supercedes
t_end, n_steps)
-
netfile
string
Name of .net file used for simulation
-
atolb
float
Absolute error tolerance for ODE’s
1e-8
rtolb
float
Relative error tolerance for ODE’s
1e-8
steady_stateb
0/1
Perform steady state check on species
concentrations
0
sparseb
0/1
Use sparse Jacobian / iterative solver (GMRES)
in CVODE
0
readFile
file
Read a .bngl or a .net file
string
Name of file to read
required
writeNET / writeSBML / writeMfile
Write current model/network in specified format
setConcentration(species,value)
Set concentration of species to value.
setParameter(parameter,value)
Set parameter to value.
saveConcentrations()
Store current species concentrations.
resetConcentratons()
Restore species concentrations to value at point of last
saveConcentrations command.
see http://bionetgen.org for BioNetGen Chapter
Simulation results with and without
equilibration
Higher initial level of
Grb2-Sos1 complex
No effect on Sos1
recruitment
More topics
• Parameter scans
• Alternate rate laws
• Labels
– Fluorescent labeling
– Carbon-fate maps
• Oligomerization
– simulate_nf New
• Compartments
• GETBONNIE
Detailed model of TLR4 signaling
• Main goals
– Represent known players and interactions in TLR4mediated signaling
– Does model qualitatively reproduce observed
tolerance in response to repeated stimulation
(“preconditioning”)?
– Investigate molecular mechanisms leading to
tolerance
G. An and J. Faeder, Math. Biosci., to appear.
Approach
• Dynamic Knowledge Representation by nonmathematician
• Relatively Detailed Representation of Components
• Mechanistic/causal emphasis
• Known source and stimulus for production of
compound (particularly inhibitors)
• Abstraction of complex molecular events if output
relatively straightforward
• Model parameters (component numbers and reaction
rates) qualitatively set based on “realistic” levels
Detailed model of TLR4 signaling
Detailed model of TLR4 signaling
Rules governing assembly of hub
complexes
Interactions between TLR4-MAL-MyD88 Complex and TRAF6
Interactions are ordered to reduce combinatorial complexity
Regulation of A20 expression
Binding of NFkB in the nucleus to A20 promoter
NFkB(Transcription~No,Activation~Yes,Location~Nucleus)+DNA(A20) \
<-> \
NFkB(Transcription~Yes!1,Activation~Yes,Location~Nucleus).DNA(A20!1)\
NFkB_DNA_A20_Bind,NFkB_DNA_A20_Unbind
Transcription of A20 gene
DNA(A20!0).NFkB(Transcription~Yes!0,Activation~Yes,Location~Nucleus) \
-> \
DNA(A20!0).NFkB(Transcription~Yes!0,Activation~Yes,Location~Nucleus) \
+ A20mRNA(Translation~On) \
A20_Transcription_Execute
Translation of A20
A20mRNA(Translation~On) -> A20mRNA(Translation~Off) + A20(TRAF6) \
A20_Translation_Execute
Degradation of A20
A20(TRAF6)->Trash(c) A20_Degrade
Generating the model
/Users/faeder/shared/Projects/BioNetGen/Perl2/BNG2.pl
BioNetGen version 2.0.47+
readFile::Reading from file TLR4_V15_PC.bngl
Read 97 parameters.
Read 31 molecule types.
Read 19 species.
Read 41 reaction rule(s).
Read 19 observable(s).
Iteration 0: 19 species
0 rxns 0.00e+00 CPU s
Iteration 1: 24 species
5 rxns 4.00e-02 CPU s
Iteration 2: 27 species 12 rxns 3.00e-02 CPU s
Iteration 3: 28 species 17 rxns 3.00e-02 CPU s
Iteration 4: 29 species 19 rxns 1.00e-02 CPU s
Iteration 5: 31 species 22 rxns 6.00e-02 CPU s
Iteration 6: 38 species 32 rxns 1.00e-01 CPU s
Iteration 7: 46 species 58 rxns 3.60e-01 CPU s
Iteration 8: 53 species 96 rxns 5.40e-01 CPU s
Iteration 9: 57 species 128 rxns 6.00e-01 CPU s
Iteration 10: 59 species 145 rxns 4.30e-01 CPU s
Iteration 11: 62 species 148 rxns 1.00e-02 CPU s
Iteration 12: 65 species 153 rxns 2.00e-02 CPU s
Iteration 13: 71 species 168 rxns 5.00e-02 CPU s
Iteration 14: 76 species 192 rxns 8.00e-02 CPU s
Iteration 15: 76 species 202 rxns 5.00e-02 CPU s
Dose-dependence of TNF-α response
and tolerance
Initial Dose Response
LPS at 10, 100, 1000, 10000
“Tolerance” Response
LPS 10000 Redose at 27h
A20 response
Initial Dose Response
LPS at 10, 100, 1000, 10000
“Tolerance” Response
LPS 10000 Redose at 27h
A20 ‘memory’ is not mechanism for tolerance
Shutting off NFkB production requires
IkB
Wild-type TNF production in
response to LPS=10,000
TNF production is not
attenuated in ‘IkB knockout’
even though A20 is produced
A20 knockout shows normal
attenuation of response
Wild-type TNF production in
response to LPS=10,000
A20 knockout