Multistate Modeling and Simulation for Regulatory Networks Zhen

Multistate Modeling and Simulation for Regulatory Networks Zhen Liu, Clifford A. Shaffer, Umme Juka Mobassera, Layne T. Watson, and Yang Cao Department of Computer Science Program in Genetics, Bioinformatics, and Computational Biology Virginia Tech

Goal: Modeling the Cell Cycle (John Tyson) G 1 l el c n o i is iv d S DNA replication M (mitosis) G 2

Regulatory Network Modeling n Model using a series of chemical reactions. The actors are proteins (“chemical species”) whose interaction rates are modeled by rate laws n Species are created, consumed, combined n Populations can rise and fall, under the control of other species n Loops and cycles n

Decomposition of Models n Modelers find it natural to divide into “bundles” of reactions.

Multistate Phosphorylation Motif n n Blocks relate to naturally occurring motifs Example: antagonistic interaction between Clb 2 and Cdh 1, with Cdc 14 as the control variable driving phosphorylation of Cdh 1 n Forms a bi-stable switch

Multistate Version n The reality is more complex, as a protein can undergo multiple levels of phosphorylation, which can affect the behavior of the larger system

Multistate Modeling n n Equations on chemical species with multiple states, related in some meaningful way Expressing as single-state equations would require dozens of reactions.

Jig. Cell Model Builder Support

Problems n n n Complications arise from the potential combinatorial explosion of states in complexes Example: Two multistate species each with 10 states could form complexes with potentially 100 states. A{i} + B{j} -> AB{i, j} This presents challenges to simulation.

Stochastic Simulation n Reaction models have often been modeled using ODEs n n Track concentrations of chemical species ODE models cannot account for stochastic effects Small numbers for some species (RNA) n Variations in inputs => Differing outputs n Simulation ensemble => Distribution n

Gillespie’s SSA (1) n n n N molecular species {S 1, …, SN}. M reaction channels {R 1, … RM}. For reaction channel Rj: Propensity function aj n State change vector vj = (v 1, j, …, v. N, j) n aj(x)dt gives probability that one Rj reaction will occur in next infinitesimal time interval given state vector x. n

Gillespie’s SSA (2) n n n Select two random numbers r 1 and r 2 Let a 0(x) be the sum for all the reaction propensities on state vector x. Time for next reaction to occur is t + t t = 1/a 0(x) log (1/r 1).

Gillespie’s SSA (3) n n n Index j for next reaction is given by smallest integer satisfying S al(x) > r 2 a 0(x). System state updated after each reaction, including populations and propensities Observations: A population-based simulation n SSA calculates propensities for reactions n

Rule-Based Modeling n n n A rule defines how a molecular particle reacts with other particles k Aopen, ? + B ---> AB, ? Subscripts describe the matching configurations for binding sites Convenient for representation Updating propensities of rules faster(? ) than updating propensities of reactions rule

Network-Free Algorithm (1) n n n (Sneddon et al. 2008) Alternative to turning rules into collections of reactions and performing SSA. Conceptually similar to SSA, but Calculate propensities for rules. n Particle based (not population based) n Keep list of particles associated with each rule n

Network-Free Algorithm (2) n Simulation loop: Calculate propensity for each rule (cheaper than SSA) n Calculate rule and time of next event n Select particles from associated list n Update the particle lists as necessary (major expense) n

Population-Based NFA (PNFA) n n (Our first contribution) Modification to NFA: (go back to) using populations for single-state species n n n Hybrid particle/population approach Attempts to cut down on the size of the lists associated with the rules Can be viewed as an optimization to NFA n at worst degrades to NFA

Full-Scale SSA (FSSSA) (1) n n (Our second contribution) Use populations even for multi-state species Should work well unless there is a small population spread across many states Can view as more direct conversion of SSA to rules (pure population-based approach)

Full-Scale SSA (FSSSA) (2) n n n For each species, store an array of populations (one for each state) Might be a sparse array Store with each rule the population count for all associated reactants

Full-Scale SSA (FSSSA) (3) n Simulation loop: Calculate propensity for each rule (cheaper than SSA) n Calculate rule and time of next event n Select a state for each reactant from the population array n Update populations of affected species (states) and population counts for associated rules (might require modifying arrays) n

Comparisons: Selection n SSA does linear search through reactions NFA, PNFA do linear search through rules, then select qualifying objects from associated reactant lists FSSSA does linear search through rules, only needs to search state lists (populations)

Comparisons: Update n n SSA updates populations of some reaction’s reactants and products NFA must create/destroy molecule objects, and update associate rule lists PNFA same, but does little work on singlestate species populations FSSSA updates sparse matrix info.

Bi-stable Switch Model n Reaction-based form: 12 species n 44 reactions n n Rule-based form: 1 single-state species, 1 multi-state n 7 rules n Non-zero populations in each state n

Simulation Times: Switch Total CPU Propensity Reactant Time Update Selection System Update Other SSA 115 72. 0 30. 6 5. 3 7. 1 NFA 341 11. 1 34. 0 286. 0 9. 9 PNFA 246 9. 9 26. 2 200. 8 9. 1 FSSSA 117 9. 2 32. 4 66. 2 9. 2

Cell Cycle Model n Reaction-based form: n n 58 species, 185 reactions Rule-based form: 17 single-state species, 6 multi-state n 64 rules n Half the states have zero population n Observation: Affecting one multi-state species affects only a smaller fraction of all the rules n

Simulation Times: Cell Cycle Total CPU Propensity Reactant Time Update Selection System Update Other SSA 171 143. 3 23. 5 1. 4 2. 8 NFA 133 36. 4 20. 4 72. 5 3. 7 PNFA 113 34. 0 17. 6 58. 6 2. 8 FSSSA 64 32. 8 18. 2 10. 5 2. 5

Simulation Quality (1)

Simulation Quality (2) n n This graph shows distribution of population for Clb 2, one of the species in the cell cycle model. The significance is that it indicates that each simulation algorithm gives approximately the same ensemble of outputs.

Complexity Analysis