Project topics and materials
Get Complete Project Material File(s) Now! »
Features of Dynamical Models in Systems Biology
We have previously detailed the crucial role that quantitative mathematical (or dynamical) models plays in modern biological reasoning. Herein, we elaborate on their required components.
The primary components of a dynamical mathematical model of a biological system corre-spond to the molecular species present in the system.The abundance of each species is assigned to a state variable of the model, the collection of which represents the system state. The system state provides a complete description of the system’s condition at any given time, via its time course, as given by the model’s dynamic behavior.
Besides state variables, models of biological systems also include parameters, which charac-terize environmental effects and interactions among system components, and whose values are fixed – meaning that the distinction between model parameters and state variables is clear-cut. Example of common parameters include association constants, maximal expression rates, and degradation rates. As a change in the value of a model parameter corresponds to a change in the environmental conditions or in the system itself, model parameters are typically held con-stant during a simulation. Varying parameter values between rounds of simulation allows one to explore system behavior under perturbations of the experimental conditions, or in altered environments.
As the main feature of a dynamical model of a biological system is its ability to describe the temporal evolution of the system components’ quantities, building a model then involves two important choices: (i) how to represent the model structure (i.e., determine its species, parameters, and molecular interactions), which is akin to defining its syntax, and (ii) how to ’interpret’, or execute, the model, which is akin to defining its semantics. These two steps are indicative of a distinction that is commonly made when dealing with systems modeling: the qualitative approach, respectively the quantitative approach. While the latter essentially relies on mathematical equations describing the performance of the system for a large set of input functions and initial states, the former requires no such formal mathematical formulation, instead relying on visual representations (e.g., diagrams) of the relationships between the system components, that is, it relies on the structural model.
For a large class of biological systems, what one usually seeks to model are the interactions between individual molecules which are only distinguishable by the class of species they belong to. Consequently, population models prove to be a well-suited, and popular, mathematical modelling framework for biological phenomena. Any population model can be described as a set of reactions operating on a set of biochemical species: this is what we will refer to as Biochemical Reaction Networks; they will represent the modeling framework of choice throughout this thesis. Below, we motivate our choice.
Network models of biological systems
While using mathematical models in the natural sciences (particularly in physics) is not a new idea1, the novelty lays in the use of such models in life sciences, most notably in biology and biomedicine. As such, quantitative mathematical models have lain outside mainstream research approaches during the last decades of the purposely reductionist qualitative era of molecular biology, which instead heavily employed qualitative models. Indeed, biologists have long used these type of models as abstractions of reality, be it in the familiar form of the ball-and-stick model of chemical structure, or more generally in the form of diagrams that illustrate a set of components and their interactions, and which play a central role in representing our knowledge of cellular processes [99].
Figure 1.1: An interaction diagram, inspired from [99]. Species A and B bind reversibly, forming a complex that inhibits the rate at which species C is transformed into species D. The graphical conventions used are that the blunt arrows denote inhibition, normal arrows denote activation, while dashed lines denote regulatory interactions, i.e. the species is not consumed by the reaction.
Such interaction maps, also called networks, formalize the complex interactions between heterogeneous entities (within and between cells), which characterize biological systems. For example, Figure 1.1 shows molecular species A and B binding reversibly, to form a complex that inhibits the rate at which species C is transformed into species D.
Generally speaking, such network-based modeling approaches help structure and formalize existing knowledge, as well as predict system behaviour, and have been widely used for describ-ing cellular regulation and metabolism [16]. Consequently, it makes sense that a quantitative modeling framework be developed using this widely-employed qualitative description of biologi-cal systems: enter Biochemical Reaction Networks (BRNs). The central argument that justifies adding quantitative/mathematical information to simple interaction maps, in order to obtain quantitative BRN models, is that simple interaction diagrams leave ambiguities with respect to the system’s behaviour: we know what components of the system interact with each other, but not necessarily how. By employing a mathematical description of the system, this uncertainty can be eliminated, at the cost of demanding a quantitative characterization of every interaction depicted in the diagram [99].
For example, in order to quantify the interaction between molecular species A and B, in Figure 1.1, a numerical description of the process must be provided, under the form of the bind-ing and unbinding reaction rate constants. For cellular processes of which only a qualitative understanding of the underlying molecular interactions is available, such a quantitative descrip-tion is not possible. However, for an important number of well-studied mechanisms, sufficient data have been collected as to allow a quantitative characterization [99], that, when coupled with the interaction diagram, can be used to formulate a mathematical model of the network’s dynamics, which typically involves the physical and chemical laws governing the mechanisms that drive the observed behaviour (such models are dubbed mechanistic). Whenever the quan-titative information regarding molecular interactions is available, the result is a collection of biochemical reactions involving the system’s species. A species is a system entity that is quan-tifiable, and whose quantity is susceptible of evolving over time. A reaction is an elementary action of the system that consists of consuming (respectively producing) a finite and integer quantity of a finite subset of species, called reactants(respectively products).
For example, the reversible binding reaction between species A and B of the interaction diagram in Figure 1.1 writes as: κ1 (1.1) −− A+B )*A.B, −−1 κ − where A.B denotes the resulting complex.
We note that, at this level of modeling, the notation A.B for the complex species is simply a name, and should not be considered as indicative of the reaction’s mechanistic details (i.e., A binds B). In this sense, one could choose any other name for the complex species: e.g., reaction κ1 1.1 could be equivalently written as A + B −−)−−* C. κ−1
The model parameters are the binding, respectively unbinding, reaction constants κ1 and κ−1. As we will see in this chapter, both the interpretation of the reaction constants and the representation of the system state will depend on the chemical kinetics assumed to govern the dynamics of the system.
Example 1.2.1
A more biologically relevant example is the Michaelis-Menten mechanism, which consists of an enzyme, denoted E, that reversibly binds a substrate, S, to form a complex, E : S. The complex then releases a product P , while preserving the enzyme. The associated reaction network writes as: κ1 κ 2 E + P (1.2)
Biochemical Reaction Networks will be the modeling framework of choice throughout this thesis. Consequently, in the next sections of this chapter we formally define both their syn-tax, as well as the two major applicable semantics: classical (deterministic) chemical kinetics, respectively stochastic chemical kinetics.
Before moving forward with the definition of BRNs, we note that a disadvantage of network models is that as the number of components and interactions in a biological system grows, it becomes increasingly difficult to maintain an intuitive understanding of the overall behaviour [99]. Thus, different levels of information abstraction can be employed when constructing network models, according to the size of the system, the available data, and the research question under consideration. In their simplest form, the interaction maps only depict the possibility of interaction between genes or proteins (Figure 1.2, top) – and can therefore cover larger systems -, while at the other end of the scale one finds detailed models based on mathematical descriptions (Figure 1.2, bottom), which provide significantly more detailed insight into the dynamics, but are usually only used to describe small, well-studied subsystems. Figure 1.2 shows a hierarchy of different, widely-employed, network models, each one abstracting information on a different level, and thus requiring different amounts of detail.
Biochemical Reaction Networks: Syntax
Models of cellular phenomena often take the form of interaction diagrams, as in Fig.1.1. For biochemical and genetic networks, interaction diagrams depict the molecular species in the system – which could be ions, small molecules, macromolecules, or molecular complexes – as nodes, and the interactions between them as arrows. The arrows can represent a range of processes: chemical binding or unbinding, reaction catalysis, or regulation of activity. The processes result in the production, inter-conversion, transport, or consumption of the species within the network.
A set of reactions constitutes a biochemical reaction network. The manner in which the biochemical species interact is referred to as the network topology, and its organization is ap-parent if one rearranges the reactions in the form of an interaction graph. For example, the interaction graph of the Michaelis-Menten mechanism of Ex.1.2.1 is shown in Fig.2.1.
In order to obtain a quantitative description, one needs to know the rates at which the reactions occur. In vivo, the rate of a reaction depends on its kinetic constant, on the con-centration of the reactants, and on physico-chemical conditions, such as temperature and pH. However, for in silico models, the physico-chemical conditions are fixed, such that rate laws can be described solely in terms of reactant concentrations and the kinetic rate constants.
All in all, a biochemical reaction network can be defined as:
Definition 2.1
(Biochemical Reaction Network) A biochemical reaction network (BRN) is a pair (S, R, α, β), such that:
(i) S = {S1, S2, . . . , Sn} is a finite set of chemical species,
(ii) R = {r1, r2, . . . , rm} is a finite set of reactions. Each reaction is a triple rj ≡ (αj, βj, κj) ∈ Nn × Nn × R≥0, that writes as:
κj + βj2S2 + . . . + βjnSn,
rj : αj1S1 + αj2S2 + . . . + αjnSn −→ βj1S1
with vectors αj and βj being commonly referred to as the consumption, respectively the production vectors of reaction rj, and κj denoting the kinetic constant of rj, (iii) α, β ∈ Rm×n denote the network’s consumption, respectively production, matrices; each element αji, respectively βji, denotes the quantity of species Si begin consumed, respectively produced, by reaction rj.
Remark 2.1. The reaction network can be written compactly in matrix-vector form, as αS −→ βS, with S the species column vector, and κ the rates column vector.
Remark 2.2. The laws of thermodynamics state that all chemical reactions are re-versible. Nevertheless, under certain environmental conditions (temperature, pressure, the availability of an enzyme), the reverse reaction proceeds at a negligible rate: nearly all of the reaction’s reactants are used to form products, which makes it very difficult, even under extreme conditions, to reverse the reaction. In this case, it is reasonable to describe the reaction as being irreversible: this is why, in Definition 2.1, we assume elementary reactions to be irreversible. Reversible reactions will simply be represented by the couple of direct (forward) and reverse (backward) reactions, which we will note, for convenience, as: kj+ αj1S1 + αj2S2 + . . . + αjnSn −−)−−* βjiS1 + βj2S2 + . . . + βjnSn. kj−
Definition 2.2
(Reactants, products, modifiers) In a reaction rj, a species Si is said to be:
• a reactant, if αji > 0
• a product, if βji > 0.
Additionally, species Si is usually considered to be a modifier in rj if it participates in the reaction both as a reactant and as a product, but is neither produced nor consumed by it: αji = βji > 0.
Definition 2.3 (Stoichiometry matrix) The stoichiometric matrix of a BRN (S, R, α, β) is formed by the stoichiometric coefficients of the reactions that constitute the network, and is defined as the m × n matrix r = β − α.
Definition 2.4 (State-change vector)
Let (S, R, α, β) be a biochemical reaction network. The jth line of the stoichiometry matrix r is called the state-change vector of reaction rj ∈ R: νj ≡ (rj1, rj 2, . . . , rjn). The state-change vector denotes the change in the molecular population caused by one occurence of reaction rj, i.e., if the system is in state x and one reaction rj occurs, the system immediately jumps to state x + νj.
Remark 2.3. We note that the process of converting a reaction network into a stoi-chiometry matrix is lossy. Otherwise said, a reaction scheme is not uniquely defined by its stoichiometry matrix, which means it is not always possible to recover the original reaction scheme from a stoichiometry matrix. For example, while the reaction systems X → ∅ and 2X → X have the same associated matrix, their dynamical behavior is clearly different.
Biochemical reaction networks represent a quantitative description of the system, which can be used to construct dynamical mathematical models. The dynamic behavior (kinetics) of a biochemical reaction network describes the changes in system states over time. A system state is described by a set of variables at a given point in time. One state should contain enough information to predict all possible future behaviours. A state can be represented in different ways, and each model defines the contents of the system state, as well as the underlying process type governing the system’s dynamics.
Generally speaking, the dynamics of a population model can be:
(i) either discrete, or continuous depending on whether the population quantity is modeled as discrete or continuous variable; discrete time is favored when the quantities of the model variables change only when specific events occur, while continuous time is favored when model variables are in constant flux;
(ii) either deterministic, if the output trajectory is fully determined by the initial state of the system, or stochastic, if starting from a given initial state, different trajectories can emerge, each trajectory having an associated probability of happening.
The notion of determinism, which is akin to model behavior reproducibility, is a foundation for much of scientific investigation. A model is called deterministic if its behavior is exactly reproducible. Although the behavior of a model depends on a specified set of conditions, no other factors influence it, so that repeated simulations under the same conditions are always in perfect agreement (i.e., they are perfect replicates).
In contrast, stochastic models allow for randomness in the model behavior, which is in-fluenced both by specified conditions and by unpredictable forces/noise. Consequently, each repetition of a stochastic simulation, under identical environmental conditions, will yield a distinct sample of system behavior.
As an illustrative example, consider the simple monomolecular reaction: κ (2.1) A−→B
The process underlying this reaction can be assumed to be either deterministic, or stochas-tic. According to both the level of abstraction and to the hypothesized nature of change in the system state in time, one further distinguishes between discrete and continuous deterministic models.
Deterministic processes with a discrete change of system state can be represented using Boolean models, which approximate the dynamics of biological networks by considering each molecule (e.g., gene or protein) in the network as either active/expressed (1) or inactive/not expressed (0). Under this representation, the questions that can asked about the system are of a qualitative nature: for example, instead of inquiring about the exact quantities of the chemical species present in the network, one is interested in their presence or absence. Despite their sim-plified underlying view of biological networks, and the fact that they introduce a coarse approx-imation by neglecting intermediate states, Boolean approaches give a meaningful insight into biological knowledge, having proven useful in analyzing system dynamics and reasoning about the stability and robustness of biological systems. For example, Boolean approaches enable the detection of singleton attractors (i.e., fixed points), where the system is stable. Morevoer, Boolean networks have been successfully applied [180] in modeling gene regulatory and signal-ing networks in a variety of biological systems ([2],[73],[118],[133],[162],[163],[166],[181]), at both cellular and population levels ([27],[116]).
However, in this thesis, we will not address discrete deterministic models such as Boolean networks. Instead, we will focus on continuous deterministic models.
For deterministic models with continuous change in time, the state of the system no longer tracks the activation state of species, but rather their concentrations. Each species in the diagram is assigned a single state variable, [Si](t), denoting its concentration at time t. The collection of values of all single state variables, {[S1](t), [S2](t), [S3](t), . . .}, at any point in time t, constitutes the state of the system. Then, to each molecular species, corresponds an ordinary differential equation (ODE) that describes how its concentration changes over time, due to interactions with other species in the network. A reaction can have one of several effects on the molecular species partaking in it: besides the obvious phenomena of species consumption and production, there are also reactions that have a modifying effect on its reactants, such as autocatalytic production, (de)phosphorylation, or (un)binding. All of these reaction are thought of as elementary: they are reactions with a single mechanistic step. A general rule for constructing the differential equation of a species S, with respect to the type of reactions it participates in, is [43]:
d[S] = synthesis − degradation − phosphorylation (2.2) dt + dephosphorylation − binding + release, etc… Equation 2.2 indicates that each reaction in which species S is consumed or structurally changed (e.g., phosphorylated, or bound to another species) contributes negatively to the evolu-tion of its concentration, while the reactions that either produce S, or revert structural changes, contribute positively to its evolution.
The rate of each reaction must be represented by a kinetic rate law, which will have one or more kinetic rate constants associated with it. In the deterministic case, these rate constants denote the speed of the reaction: how frequently the reaction is expected to occur. As we will see in the next section, biochemical reactions are assumed to operate under the law of mass action, which states that the rate of a chemical reaction is directly proportional to the product of the activities or concentrations of the reactants.
For example, under the law of mass action, the changes in concentration for species A and B in Example 2.1, in the time interval dt, are given by the following system of ODEs:
(d[A] = −κ[A] dt
d[B] = κ[A] dt
which states that species B is produced at the same rate at which species A is consumed, and that that rate is proportional to both the concentration of the reaction’s reactant A, and to its kinetic rate constant κ.
The third modeling choice assumes that the underlying process is non-deterministic, i.e., that it exhibits a random component. Then, one can employ stochastic models, in which the state system is given by a vector of species’ molecule count (xS1 (t), xS2 (t), . . .) ∈ N≥0, instead of concentrations, and whose dynamics is governed by the probability of a reaction occurring in a small time interval (t, t + δt]. In this case, the rate constant of a reaction rj is no longer interpreted as the reaction’s speed, but rather as an indicator of the probability that a tu-ple of molecules corresponding to rj’s reactant species will react according to rj in the next infinitesimal time dt.
Then, the rate of reaction A −→ B denotes the probability of a molecule A transforming into a molecule B in a time interval dt: P(xa(t + dt) = xa(t) − 1, xb(t + dt) = xb(t) + 1) = κxa(t), with (xa(t), xb(t)) denoting the number of molecules of type A, respectively of type B, present in the system at time t.
As the latter two frameworks constitute the modelling approaches used in this manuscript (i.e., we do not tackle Boolean models), we next formally elaborate on the concepts of stochastic and deterministic semantics of BRNs.
Classical chemical kinetics
Conventional chemical kinetics operate under a continuous deterministic modelling frame-work, which involves no randomness in the development of states. That is, if a deterministic system is known at one time, then it is known at all subsequent times. It serves as the traditional way of representing the systems dynamic of complex biological systems that involve the interac-tion of many components: starting in the late 1970s, researchers began modeling cell physiology primarily using the continuous deterministic ODE approach, and creating increasingly detailed models over the next three decades.
Under this framework, the biochemical reaction network is modeled as a continuous system with continuous time dynamics. For a given BRN (S, R, α, β), the system state is represented by a n-vector (x1(t), x2(t), . . . , xn(t)) of continuous variables that keep track of reactant con-centrations, i.e., xi(t) denotes the concentration of species Si in the system at time t. The chemical reactions/interactions are also represented by continuous processes, i.e., it is assumed that reactions occur continuously and simultaneously.
Each process rj has an associated deterministic rate constants kj – whose value is identical to that of the dimensionless constant κj introduced in Definition 2.1. It gives a measure of how frequently each type of reaction is expected to occur. The velocity of each reaction is specified using a rate equation. Deterministic reaction rates are described under the following assumptions:
(i) stochastic fluctuations in the system are negligible;
(ii) molecules are considered to be point-like entities;
(iii) the reaction volume border, as well as the resulting frontier effects, are neglected;
(iv) spatial homogeneity: the reaction volume is well stirred, i.e., the reactants are equally distributed throughout the volume, meaning that the rate of each reaction is independent of the reactant position in space; the reaction rate can then be unambiguously referred to in the volume;
(v) the continuum hypothesis: each species is present in the system in a large amount – molec-ular abundance can be then described in terms of their continuously varying concentration, as opposed to an integer-valued molecule count.
Then, in a fixed volume, and under the negligible stochastic fluctuation, spatial homo-geneity and continuum hypothesis, the reaction rate equations typically assume mass-action kinetics – or an enzyme kinetic law based on mass-action, such as the Michaelis-Menten or Hill kinetics. The law of mass-action states that the rate of a chemical reaction is proportional to the product of the reactants’ concentration.
Under these conditions, the dynamic evolution of the system state (represented as an ar-ray of species’ concentrations) can be described mathematically by a set of coupled ordinary differential equations, called the reaction-rate equation (RRE):
Definition 3.1.1 (Continuous deterministic model (RRE))
Let (S, R, α, β) be a biochemical reaction system and x0 = (x1, x2, . . . , xn) ∈ Rn≥0 an initial state of the system. Then, the continuous deterministic model is the solution of the set of n coupled differential equations given by1: dx(t) = rT f(x(t)), (t ∈ R≥0) (3.1)
Table of contents :
Introduction
Computer Science and Systems Biology
Challenges
Outline and Contributions of the Thesis
I Biochemical Reaction Networks: A Review
1 Context and Motivation
1.1 Features of Dynamical Models in Systems Biology
1.2 Network models of biological systems
2 Biochemical Reaction Networks: Syntax
3 Biochemical Reaction Networks: Semantics
3.1 Classical chemical kinetics
3.2 Stochastic chemical kinetics
3.3 Stochastic vs deterministic models
4 Executable Biology: Computational models of Biochemical Reaction Networks
4.1 Petri Nets
4.1.1 Motivation
4.1.2 Qualitative Petri Nets
4.1.3 Quantitative Petri nets
4.2 Rule-based modeling and the Kappa language
4.2.1 Rule-based modeling
4.2.2 The Kappa language: Syntax and Operational Semantics
4.2.3 The Kappa language: Stochastic Semantics
II Model reduction
Motivation
5 Prototyping genetic circuits using rule-based models
5.1 Motivation
5.2 A Kappa model of the -phage decision circuit
5.2.1 Example: Kappa model of CI and CII production
5.3 Model Approximation using Michaelis-Menten-like reaction schemes
5.3.1 Validity of the Michaelis-Menten enzymatic reduction in stochastic models
5.3.2 Generalized competitive enzymatic reduction
5.3.3 Operator site reduction
5.3.4 Fast dimerization reduction
5.3.5 Top-level reduction algorithm
5.4 Results and Discussion
5.4.1 Results: reduction of the -phage decision circuit
5.4.2 Conclusion and future work
6 Tropical Abstraction of Biochemical Reaction networks with guarantees
6.1 Introduction
6.2 Definitions and Motivating Examples
6.2.1 General Setting and Definitions
6.2.2 Motivating example: Michaelis-Menten
6.2.3 Motivating example: A DNA model
6.3 Model reduction using conservative numerical approximations
6.4 Error estimates of tropicalized systems using conservative numerical approximations
6.4.1 Tyson’s Cell Cycle Model
6.5 Comparison with existing methods
6.6 Conclusion and outlook
III Modelling storage and growth
7 Effects of cellular resource storage on growth
7.1 Introduction
7.2 Review of the Weisse cell model
7.2.1 Overview
7.2.2 Nutrient Uptake and Metabolism
7.2.3 Transcription
7.2.4 Competitive Binding
7.2.5 Translation
7.2.6 Growth and cellular mass
7.2.7 Model parameters
7.3 A scaling procedure for modelling nutrient storage in the cell
7.3.1 Motivation
7.3.2 Definition
7.3.3 Examples
7.3.4 Storage capacity in the Weisse model
7.4 The Growth Rate is an emergent property of the model
7.5 The Single Cell model
7.5.1 The effects of metabolite storage on exponential growth rate
7.5.2 Metabolite storage and adaptation to environmental fluctuations
7.5.3 Environmental dynamics and resource storage strategies
7.6 The Ecological perspective
7.6.1 Motivation
7.6.2 Experiment setup: a mixed population model in a fluctuating environment
7.6.3 Results
7.7 Conclusion
8 Abstractions for cellular growth: towards a new Petri net semantics
8.1 Motivation
8.2 Split-Burst: towards a piecewise-synchronous execution of Biochemical Reaction Networks
8.2.1 Max-parallel execution semantics of Petri Nets
8.2.2 Piecewise-synchronous execution semantics
8.3 Linear Reaction Networks
8.3.1 Growth rate as an emergent property
8.3.2 Synchronous execution
8.4 Possible Applications
8.4.1 Approximation of system dynamics
8.4.2 Synchronous Balanced Analysis
8.5 Conclusion and future directions
Conclusion
9 Conclusion and future directions
9.1 Contributions
9.2 Future works
Appendices
A Probability Theory
B A DNA model: equation for the derivative of the lower bound on the concentration of x2
Bibliography
