Документ взят из кэша поисковой машины. Адрес оригинального документа : http://www.iki.rssi.ru/microsoft/ppt/mura.pdf Дата изменения: Tue Jun 30 16:25:59 2009 Дата индексирования: Tue Oct 2 06:13:11 2012 Кодировка:

CoSBiLab: a software framework to support incremental modeling
Ivan Mura
mura@cosbi.eu

Centre for Computational and Systems Biology, Trento ­ Italy


CoSBi in numbers

Staff: total 31
Sr. Res Jr. Res Developers 2 Res Ph.D. St Admin
Italy Japan

Nationality
Hungary Rumania UK Vietnam
Math

Area of research
Computer Science Life Sciences Phys

USA 1
1

Turkey 11

8

3

1 2

2

1

5 16
21

2

3

9
7

IKI meeting, 10-12 June 2009


CoSBiLab


What we do at CoSBi In-silico lab to support modeling through refinement
Modeling tools at work




specification



simulation
analysis/validation refinement





IKI meeting, 10-12 June 2009


Specification

IKI meeting, 10-12 June 2009


Fuzzy pairing in biology
Key-Lock model 1 key opens 1 lock

Computer science
· interactions occur as communications on a channel

exact matching name/co-name

Master-Key model 1 key opens a set of locks

Observed biological phenomena
· protein interactions · DNA binding

fuzzy matching, variable affinity

IKI meeting, 10-12 June 2009


BlenX
Biological transformations as interactions of bio-processes equipped with affine binders
· bio-processes inspired to pi-calculus · typed binders, affinity defined between type pairs

One box, intuitively
Typed interaction sites Kernel of the box: two main processes running in parallel Message_Handler: manages interaction protocols between boxes Message_Handler || Interface_Handler

Interface_Handler: manages the modifications of the typed interaction sites

IKI meeting, 10-12 June 2009


Language flavor
Binders

let P : pproc = x?(m).y!(m).z!(m).nil; let Box : bproc = #(x,A),#(y,B),#(z,C)[ P ];

Internal process

Stochastic · each interaction happens in a random time · the average time of interaction is determined by the affinity
IKI meeting, 10-12 June 2009


The BlenX language: elementary reactions
Interaction capabilities

Domain Species
P Internal behaviour

Monomolecular reaction
IEP IEP
in

Dom1:IEP_A change(Dom1, IEP_IN).P IEP change

Dom1:IEP_IN
P IEP

Complexation
Cyc B Cdk Cyc B

Dom1:CycB_A Dom2:Cdk_A P CycB Q Cdk complex

Dom1:CycB_A Dom2:Cdk_A

P CycB

Q Cdk

Cdk

IKI meeting, 10-12 June 2009


The BlenX language: abstract reactions
Interaction capabilities

Domain Species
P Internal behaviour

Enzymatic activation
IEP

Dom1:IEP_A
Cdc20

Dom2:Cdc20_IN Dom2?().Q Cdc20 inter

Dom1:IEP_A
P IEP

Dom2:Cdc20_A Q Cdc20

Cdc20

in

Dom1!().P IEP

Regulated degradation
Cdh1

Dom1:Cdh1_A Dom2:CycB_A

Dom1:Cdh1_A inter P Cdh1

Dom2:CycB_A die CycB

Cdk Cyc B

Dom1!().P Cdh1

Dom2?().die CycB

IKI meeting, 10-12 June 2009


Dynamic creation of entities


Biological compounds have sites of interaction
multiple sites can be present in the same entity bindings occur reversibly between 2 affine sites complexes of biological components can assembly without a precise order and can result in different topological structures example: protein C has 2 sites, both affine to 2 sites of protein W


W

C

C

W

C

How many structures can form?



It may be cumbersome or even impossible to specify such a behaviors in many modeling formalisms

IKI meeting, 10-12 June 2009


Simulation

IKI meeting, 10-12 June 2009


Intrinsic discreteness


The truly molecular nature of biological interaction was considered hardly tractable


tracking single molecules state, location and movement is indeed quite heavy from a computational point of view



In 1976, D. T. Gillespie
proved that the evolution of a well-stirred biochemical system can be represented as a Markov process provided a very simple and extremely efficient simulation algorithm for computing realizations of such process




Gillespie's algorithm (SSA) paved the way for a number of discrete modeling approaches

IKI meeting, 10-12 June 2009


Direct method (1976)


Model


N species, M reactions X(t) represent the state of the system at time t, t0 (a stochastic process)



Given X(t)=x , the probability that the next reaction happens in the infinitesimal time interval [t+t,t+t+dt] and is a reaction of type j is aj(x) exp(-a0(x) t)
the time t to the next reaction is an exponential random variable of mean 1/a0(x) the probability that next reaction is of type j is aj(x) /a0(x)




At each simulation step, 2 uniform r.n. u and v are drawn
t is chosen to be ln(u-1)/a0(x) j is chosen as the smallest integer satisfying



i 1

j

ai ( x ) v a0 ( x )

IKI meeting, 10-12 June 2009


Reformulations of the method


First reaction method (1976)
at each simulation step, draw m uniform r.n. and compute t1, t2,..., tm, the putative time of all reactions choose t as the min(t1, t2,..., tm) choose j as the index of the minimum above




Next reaction (2000)
same as the above one, but the putative times are saved in an indexed binary tree so that the minimum is always at the top a dependency graph is used to keep track of coupling among reactions to determine when putative times in the tree have to be resampled




Modified direct method (2004)


a pre-run to determine a suitable order of reactions to minimize cost of step 2) self-adaptive version of the one above, no pre-run



Sorting direct method (2006)


IKI meeting, 10-12 June 2009


Rare events


Events that lead to abnormal states are rare but have devastating effects


rare epigenetic modifications play crucial roles in cancer development failure of DNA repair mechanisms occurs randomly with a very low probability per replication





SSA computational requirements for the analysis of rare events may be substantial Example




the spontaneous switch from the lysogenic to the lytic state in phage -infected E. coli is experimentally estimated to be in the order of 10-7 per cell per generation SSA would generate sample trajectories of this rare event once every 107 runs


it would require more than 1011 simulation runs to generate an estimated probability within 1% relative half-width of a confidence interval (at a 95% confidence level)

IKI meeting, 10-12 June 2009


Biasing SSA trajectories: wSSA


Trajectories generated by the SSA are determined by the reaction propensities aj(x), j=1,2,...,M


vector field in the phase portrait
rare events





Rare event most trajectories do not reach it Biasing propensities (I.S. technique) can increase the likelihood of finding it


more samples in lesser simulation runs increased precision with the same computation time bias the propensities remove the bias from results
initial state



There is a simple way to


Kuwahara, Mura. An efficient and exact stochastic simulation method to analyze rare events in biochemical systems, JCP, 129(16), 2008
IKI meeting, 10-12 June 2009


Computational saving


Example: simple synthesis/degradation system


A, with rate k1=1 A , with rate k2=0.025 analytical solution, steady state distribution of A is Poisson, average A=40



We estimate the probability that, starting at t=0 with A=40, the systems reaches a state with A={65,70,75,80} within the interval [0,100]
2.0 1.8 1.00E+16 1.00E+14 1.00E+12 1.00E+10 1.00E+08 1.00E+06 65 70 75 80

1.6
1.4 1.2 1.0 0.8 0.6 0.4 0.2 0.0

ratio SSA/wSSA to achieve a fixed accuracy

ratio wSSA/SSA for a fixed number of runs

1.00E+04 1.00E+02 1.00E+00

65 80



When the rare event is defined as "reaching the state when A=80 within *0,100+", wSSA can compute an estimate that is 99.9999% accurate in 5.6104 secs, whereas SSA would require 2.3 105 years

IKI meeting, 10-12 June 2009


Analysis/validation

IKI meeting, 10-12 June 2009


Oscillatory behaviors


Many biological systems achieve equilibrium conditions that are not commonly found in artificial systems


living systems keep oscillating

6000 5000 4000 3000 2000 1000 0 0

Stochastic solution (1 run)
Predators Preys



Many systems have transient oscillation that stop abruptly


dead

1

2

3

4

5

6

7

8

9 10



This poses issues in
defining adequate measures that can characterize cyclic system behavior comparing variants of models of systems


IKI meeting, 10-12 June 2009


Fourier analysis
Analysis of stochastic simulations for oscillatory phenomena


identification of oscillation determination of period from noisy traces

A time to frequency transformation can facilitate the job
sin(t/10) + sin (t/20) + w.n. in time sin(t/10) + sin (t/20) + w.n. in freq.

FT

IKI meeting, 10-12 June 2009


Comparison between models
Wild Type and cln1 cln2 cln3 sic1 S. cerevisiae mutant

1.0 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0.0

spectra obtained from multiple stochastic simulation traces
Mass CycB

0
1. 0 0. 9 0. 8 0. 7 0. 6 0. 5 0. 4 0. 3 0. 2 0. 1 0. 0 0

100

200

300

400

500

600

700

800

900

1,000

Mass

CycB

100

200

300

400

500

600

700

800

900

1, 000

IKI meeting, 10-12 June 2009


Comparison of models and experiments
NF-kB oscillations
Model 1

Experimental

Model 2

r = 0.0625 r = 0.0352

IKI meeting, 10-12 June 2009


Refinement

IKI meeting, 10-12 June 2009


Including new experimental data
A model gets refined by inclusion of new experimental evidence
· ·

additional species detailed reactions

We are working on a tool for estimating the kinetic parameters of models from wet-lab data samples
Time series of reactant concentrations · Rate constants · Error on rate constant estimates · Variance on rate constant estimates

Experimental data

KInfer
Structure of biochemical network Chemical reactions of the system

IKI meeting, 10-12 June 2009


ML approach
The values of the model's parameters have to be the most likely values giving the observed time course of the concentration.

The Maximum likelihood (ML) approach can be used to achieve this goal. The main steps are:
to build a suitable expression of the joint transitional density for expressing the probability density function of the observed outcomes in terms of measured system variables and parameters;
and to optimize this function to determine unknown parameters.

Lecca, Palmisano, Priami, Sanguinetti, A new probabilistic generative model of parameter inference in biochemical networks, ACM Symposium on Applied Computing 2009

IKI meeting, 10-12 June 2009


Kikuchi model, benchmarking rate inference
Kikuchi et. al, "Dynamic modeling of genetic networks using genetic algorithm and S -system", BIOINFORMATICS, 19(5) 2003

A gene regulation network with complex feedbacks Our results
·

simulations with inferred rates

IKI meeting, 10-12 June 2009