Mesoscopic and Stochastic Phenomena and the Lac operon

1. Mesoscopic and Stochastic Phenomena and the Lac operon Ádám Halász joint work with Agung Julius, Vijay Kumar, George Pappas Adam Halasz

2. Outline • Lactose induction in E. coli, an example of bistability • Stochastic phenomena in reaction networks • Mesoscopic effects in the lac operon • Outlook Adam Halasz

3. Lac system: Biological phenomenology lac Z lac Y lac A lac I mRNA mRNA b galactosidase repressor permease internal lactose allolactose external lactose Adam Halasz

4. mRNA M Permease P β-galactosidase B External lactose Allolactose A Lactose L Lac system: ODE model • Network of 5 substances • Example of positive feedback in a genetic network discovered in the 50’s • Described by differential equations which are built from chemical rate laws • Some time delays and time scale separations ignored and/or idealized Adam Halasz

5. Lac system: ODE model mRNA M • Network of 5 substances • Example of positive feedback in a genetic network discovered in the 50’s • Described by differential equations which are built from chemical rate laws • Some time delays and time scale separations ignored and/or idealized Permease P β-galactosidase B External TMG Te TMG T Adam Halasz

6. Pout Pin Bequilibrium Texternal Lac system: ODE model • Because of the positive feedback, the system has an S-shaped steady state structure • That is to say, for some values of the external inducer concentration (Te), there are two possible stable steady states P B Te Adam Halasz

7. Lac system: Bistability • Switching and memory • Need to clear T2 in order to switch up B B t Te t Te T1 T2 Adam Halasz

8. Lac system: Bistability • Switching property is robust • Model parameters perturbed by 5% Adam Halasz

9. Lac system – first lesson • Network of <5 species involved in reactions • Reactions decomposed* into mass action laws • Can be implemented using simple stochastic transition rules eg: “upon colliding with a B, A becomes C with probability x” • Recipe for synthesis of switch with hysteresis • Compose motifs to build logical functions,… • Design the network and transition rules • Map the rules to individual stochastic programs Adam Halasz

10. Lac system: ODE is not enough Adam Halasz

11. Stochastic versus deterministic • Substances are represented by finite numbers of molecules • Rate laws reflect the probability of individual molecular transitions • The abstraction/limit process is not always trivial Next transition time distribution:: t2 N2=N1-2 t1 N1 Adam Halasz

12. A + B AB A + B AB Chemical reactions are random events B B A A Adam Halasz

13. Stochastic reaction kinetics • Quantities are measured as #molecules instead of concentration. • Reaction rates are seen as rates of Poisson processes. k A + B  AB Rate of Poisson process Adam Halasz

14. Stochastic reaction kinetics A AB time reaction reaction reaction time Adam Halasz

15. k k 1 2 Multiple reactions • Multiple reactions are seen as concurrent Poisson processes. • Gillespie simulation algorithm: determine which reaction happens first. A + B AB Rate 1 Rate 2 Adam Halasz

16. Multiple reactions A AB time reaction 1 reaction 2 reaction 1 time Adam Halasz

17. When the number of molecules per cell is small*, the respective substance has to be treated as an integer variable The probabilistic transition rules can be implemented in standard ways Gillespie method: instead of calculating time derivatives, we calculate the time of the next transition Many other sophisticated methods exist. As an empirical rule, the higher the number of molecules, the closer the simulation is to the continuous version Gillespie simulations Adam Halasz

18. Gillespie automata • The state of the system is given by the number of copies of each molecular species • Transitionsconsist of copy number changes corresponding to elementary reactions • The distribution of the next transition timeis Poisson, e-kt where k is the propensity • A Gillespie automatonis a mathematical concept [a continuous time Markov chain] • Plays the same role differential equations have in the continuum description Adam Halasz

19. mRNA concentration Increase E # mRNA molecules Time (min) External TMG concentration Lac system: spontaneous transitions both ways Adam Halasz

20. Aggregate simulations Mixed Gillespie/ODE 1000 cells Adam Halasz

21. Lo Hi Two state Markov chain model • We create a simplified model, a continuous time Markov chain with two discrete states, high state and low state • The transition rates depend on the external concentration of TMG (Julius, Halasz, Kumar, Pappas, CDC06, ACC07) Adam Halasz

22. Transition rates Identified transition rates (Monte Carlo) Adam Halasz

23. Two state Markov chain model Average of a colony with 100 cells # mRNA molecules E[M ] t Time (min) Adam Halasz

24. Lac system – second lesson • Bad news: underlying stochasticity can drastically modify the ODE prediction • The price paid for an ‘uncontrolled’ approximation • Good news: the ODE abstraction can also be performed rigorously • Gillespie automaton is an exact abstraction • For Gillespie  ODE, all ‘molecule’ numbers must be large • Stochastic effects retained at the macro-discrete level • Effects are reproducible and quantifiable • Further abstractions of stochastic effects are possible • Lac example: can quantify the spontaneous transitions: • Choose an implementation where they are kept at a low rate • Implement control strategies that use the two-state model Adam Halasz

25. Beyond Gillespie • Gillespie method is ‘exact’ – produces exact realisations of the stochastic process • Main problem is computational cost • Larger molecule numbers • Rare transitions • Several approaches to circumvent ‘exact’ simulations Adam Halasz

26. Towards the continuum limit t – leaping: lump together several transitions update molecule numbers at fixed times: A AB time r2 r1 r2 r1 r1 r2 r1 D D D D time Adam Halasz

27. Towards the continuum limit • Error introduced by t – leapingisdue to variation of the propensities over the time interval (may lead to negative particle numbers!) • Acceptable if the expected relative change of each particle number over Δ is small (e.g. if the number of particles is large) • If the number of transitions per interval is also large, the variation can be described as a continuous random number  stochastic differential equations • Finally, is the variance of the change per interval can also be neglected, the simulation is equivalent to an Euler scheme for an ODE. Adam Halasz

28. Other limiting cases • If the number of all possible configurations is relatively small, probabilities for each state can be calculated directly, by calculating all possible transition rates, (finite state projection) or using the master equation (Hespanha, Khammash,..) • In some situations (eg. signaling cascades) there is a combinatorial explosion of species, where agent-based simulations are useful (Los Alamos group, Kholodenko) Adam Halasz

29. Transitions in the lac system • We used tau-leaping for our simulations • The high state can be simulated using SDEs or ‘semiclassical’ methods • The lower state can be studied using finite state projection # mRNA molecules Time (min) Adam Halasz

30. Summary • Mesoscopic effects in biological reaction networks are due to small numbers of molecules in individual cells • They may affect the system dramatically, somewhat, or not at all • These effects can be described mathematically and incorporated in our modeling efforts • Several sophisticated methods exist; it is important to use an approximation that is appropriate, both in terms of correctness and in terms of efficiency Adam Halasz

31. Thanks: Agung Julius, , George Pappas, Vijay Kumar, Harvey Rubin DARPA, NIH, NSF, Penn Genomics Institute Adam Halasz