Coloured Petri nets are imperative for studying bigger biological models, particularly, those which expose repetition of components. Such models can be easily scaled by minor changes of a few parameters. Similarly, studying certain biological phenomena necessitates the existence of discrete and continuous variables as well as continuous and stochastic processes in one and the same model. Thus, hybrid modelling and simulation is indispensable to deal with these challenges. In this paper we introduce a generic Petri net class, Coloured Generalised Hybrid Petri Nets (GHPN C ) by combining coloured Petri nets and Generalised Hybrid Petri Nets (GHPN ), which integrates discrete and continuous places as well as stochastic and deterministic transitions on the coloured level. Moreover, we present a case study which illustrates a possible and typical application where such a Petri net class is highly required.
Petri nets have been widely used for modelling and analysis of biological systems, see e.g., [2,10,18,21,31,33,34,35,36,39,40]. Their intuitive graphical representation makes them easily approachable by biologists. Furthermore, Petri nets possess well-established mathematical notations which may render them as the future de facto standard for modelling biological systems compared to other graphical languages currently in use. Nevertheless, with the rapid progress of systems biology, standard place/transition Petri nets become insufficient to accommodate the needs of systems biologists to study larger models. Therefore, many different extensions have been adapted for their potential contribution to systems biology. Among such promising extensions are hybrid Petri nets (HPN ) [1,6] and coloured Petri nets (PN C ) [23].
Hybrid Petri nets [1,6] have been increasingly motivated for their contribution to systems biology [17,18,35,36]. On the one hand (biological) systems may occur at different time scales: slow and fast [14,25]. Using one modelling paradigm (i.e., discrete or continuous paradigm alone) tends to be inefficient for studying multitimescale biological models [18]. On the other hand the integration of discrete and continuous variables as well as deterministic and stochastic processes in the same model is necessary in certain application scenarios to implement particular model semantics [18,36]. For instance, in [21] discrete places and immediate transitions were employed to implement the part related to cell division of the eukaryotic cell cycle while continuous and stochastic transitions were used to represent and quantitatively simulate biological reactions. Additionally, the HPN formalism provides a different approach to control the accuracy and the speed of biological models during simulation [18]. Thus, they provide the modeller with a tool to make a trade-off between the simulation's efficiency and the result's accuracy.
Another, yet powerful extension of standard Petri nets is coloured Petri nets (PN C ) [23,31,40], where a group of similar components are represented by one component, each of which is defined as and thus distinguished by a colour [27]. PN C are very useful for modelling larger biological systems where low-level Petri nets do not scale well, while a (biological) system modelled via coloured Petri nets can easily be scaled with minor modifications of certain coloured variables. Furthermore, modelling of biological systems is moving from single to multiple scales (multi-scale modelling), which introduces a series of challenges such as repetition of components, communication between components, organisation of components, and pattern formation of components [27]. All these challenges potentially could be tackled by coloured Petri nets rather than standard Petri nets.
Nevertheless, with the rapid change of type and size of models of biological systems which have to be analysed, a Petri net class that integrates all the features of hybrid Petri nets and coloured Petri nets is highly required. The high-level representation of coloured Petri nets can be used to model larger biochemical networks and stochastic and continuous components can be simultaneously used to facilitate the efficient simulation of multi-timescale models.
Thus, in this paper we introduce a class of Petri nets, Coloured Generalised Hybrid Petri Nets (GHPN C ) by combining all features of Generalised Hybrid Petri Nets (GHPN ) [18] and the merits of Coloured Petri Nets (PN C ) [27]. The new class supports a rich set of transition types as well as discrete and continuous places. Moreover, it permits the full interplay between stochastic and deterministic processes at the coloured level. All the feature of GHPN C are implemented in Snoopy [15] -a unifying Petri net tool which is available free of charge for academic use.
The paper is structured as follows: first we briefly discuss the related work of hybrid Petri nets and coloured Petri nets followed by a brief background of each of those net types. Next, we present the formal definition of the new Petri net class followed by its semantic as well as an outline of the simulation algorithm used to produce the dynamic of GHPN C . Afterwards, we present one possible and typical application of GHPN C in the context of systems biology, the repressilator model, which is an engineered synthetic system encoded on a plasmid. We conclude by summarizing possible extensions of Coloured Generalised Hybrid Petri Nets.
In this section, we will briefly describe related work on hybrid and coloured Petri nets for systems biology.
Hybrid Petri nets have been introduced in [1] to deal with situations where discrete and continuous entities exist in the same model. The motivation given in early publications was in modelling technical systems whereby discrete places are used to model the states of switch-like components while continuous places are used to represent their fluid counterparts. The pioneer work of introducing hybrid Petri nets to systems biology was done by Matsuno et al.,in [36]. They had noticed that using equal inflow and outflow in the fluid part of hybrid Petri nets is not natural to use them to model biological processes. Thus, they introduced Hybrid Functional Petri Nets (HFPN) [36,35]. Many successful models have been implemented using HFPN (e.g., see [8,34]).
In [18], we have used Hybrid Petri nets in a different way whereby stochastic transitions are used to model slow biological processes, while fast processes are modelled via continuous transitions. In [21], this approach is used to model the eukaryotic cell cycle.
The early applications of coloured Petri nets for systems biology were limited, which, to our knowledge, can almost only be seen in [2,3,5,11,26,38,40,41,42]. These studies are rather small and usually resort to Design/CPN [4] or its successor CPN Tools [23]. However these tools were not specifically designed with the requirements of systems biology in mind. Thus they are not suitable in many aspects, e.g., they do not directly support stochastic or continuous modelling, nor stochastic or deterministic simulation.
Since coloured Petri nets were introduced to our Petri net modelling tool, Snoopy, we have extensively explored the application of coloured Petri nets for (multiscale) systems biology. For example, in [29], we used coloured stochastic Petri nets to model and analyse stochastic membrane systems, where each compartment is encoded as a colour. In [30], we described the multiscale modelling of coupled Ca 2+ channels using coloured stochastic Petri nets by considering two levels: Ca 2+ release sites and Ca 2+ channels. In [10], Coloured stochastic and coloured continuous Petri nets are used for multiscale modelling and analysis of Planar Cell Polarity in the Drosophila wing, and the built model consists of more than 800 cells. In [37], a case study of phase variation in bacterial colony growth is given, in which cells are distributed on a two-dimensional grid represented by both Cartesian and polar coordinate systems.
In [18] we have introduced a special class of Hybrid Petri nets called Generalised Hybrid Petri Nets (GHPN ). The main two objectives of GHPN are: to provide systems biologists with a convenient and flexible graphical tool to model and simulate biological processes with different time scales and to facilitate the process of constructing hybrid models where stochastic and deterministic (i.e, continuous) processes are tightly coupled. To fulfil such objectives, GHPN contain a rich set of element types: places (discrete, and continuous), transitions (stochastic, immediate, deterministic, stochastic, and continuous), and arcs (standard, read, inhibitor, equal, and reset) [18,17]. Figure 1 depicts the graphical representation of the different element types of GHPN .
Discrete places (single line circle) hold non-negative integer numbers which may represent the number of molecules of a given species (tokens in Petri net notions). On the other hand, continuous places -which are represented by the shaded line circle -hold non-negative real numbers which represent the concentration of a certain species. Please note that, except when otherwise mentioned, the number which a place p i holds, also called its marking, is referred to by m(p i ).
Furthermore, GHPN offer five transition types: stochastic, immediate, deterministically delayed, scheduled, and continuous transitions [18]. Stochastic transitions, which are drawn in Snoopy as a single line square, fire randomly with an exponentially distributed random delay. The user can specify a set of firing rate functions that determine the random firing delay. The transitions' pre-places can be used to define the firing rate functions of stochastic transitions. Immediate transitions (black bar) fire with zero delay, and have always highest priority in the case of conflicts with other transitions. They may carry weights (which can also be defined by state-dependent functions) that specify the relative firing frequency in the case of conflicts between immediate transitions. Deterministically delayed transitions (represented as black squares) fire after a specified constant time delay. Scheduled transitions (grey squares) fire at user-specified absolute time points. Continuous transitions (shaded line square) fire continuously in the same way as in continuous Petri nets. Their semantics are governed by a set of ordinary differential equations (ODEs) which define the changes in the transitions' pre-and post-places. More details about the biochemical interpretation of deterministically delayed, scheduled, and immediate transitions can be found in [16].
The connection between those two types of nodes (places and transitions) takes place using a set of different arcs. GHPN offer six types of arcs: standard, inhibitor, read, equal, reset, and modifier arcs. Standard arcs connect transitions with places or vice versa. They can be discrete, i.e., carry non-negative integervalued weights (stoichiometry in the biochemical context), or continuous, i.e., carry non-negative real-valued weights. In addition to their influence on the enabling of transitions, they also affect the place marking when a transition fires by adding (removing) tokens from the transition's post-places (pre-places). Extended arcs such as inhibitor, read, equal, reset, and modifier arcs can only be used to connect places to transitions, and not vice versa. A transition connected with an inhibitor arc is enabled (with respect to the corresponding pre-place) if the marking of the pre-place is less than the arc weight. In contrast, a transition connected with a read arc is enabled if the marking of the pre-place is greater than or equal to the arc weight. Similarly, a transition connected using an equal arc is enabled if the marking of the pre-place is equal to the arc weight.
The other two remaining arcs do not affect the enabling of transitions. A reset arc is used to reset a place marking to zero when the corresponding transition fires. Modifier arcs permit one to include any place in the transitions' rate functions and simultaneously preserve the net structure restriction. Besides, the markings of places connected using read, inhibitor, equal, or modifier arcs does not change when the corresponding transition fires.
Coloured Petri nets [22] consist, as standard Petri nets, of places, transitions and arcs that connect places and transitions. Additionally, a coloured Petri net is characterised by a set of finite colour sets. Each place gets assigned a colour set and may contain distinguishable tokens; each token has a colour of this colour set.
As there can be several tokens of the same colour on a given place, the tokens on a place define a multiset over the place's colour set.
Each transition gets a guard, which is a Boolean expression over variables, constants, functions, etc. The guard of a transition must be evaluated to true for the enabling of the transition. The trivial guard "true" is usually not explicitly given.
Each arc gets assigned an expression, rather than an integer number in standard Petri nets; the result type of this expression is a multiset over the colour set of the connected place.
In this section we present the definition of coloured generalised hybrid Petri nets. GHPN C reuse all elements of GHPN , however on the coloured level. We start off by defining GHPN C . Afterwards, we present the semantics of GHPN C that governs the execution of models constructed by it. [27,18], where:
-P = P disc ∪ P cont whereby P disc is the set of discrete places and P cont is the set of continuous places.
1. T stoch is the set of stochastic transitions, which fire stochastically after an exponentially distributed waiting time. 2. T im is the set of immediate transitions, which fire with waiting time zero; they have higher priority compared with other transitions. 3. T timed is the set of deterministically delayed transitions, which fire after a deterministic time delay. 4. T scheduled is the set of scheduled transitions, which fire at predefined time points. 5. T cont is the set of continuous transitions, which fire continuously over time. where
is a finite, non-empty set of colour sets. -C : P → is a colour function that assigns to each place p ∈ P a colour set C(p) ∈ .
-F : A → EXP is an arc function that assigns to each arc a ∈ A an arc expression of a multiset type C(p) M S , where p is the place connected to the arc a. -V is a set of functions V = {g, d, w, f } where :
1. g : T stoch → H s is a function which assigns a stochastic hazard function h st to each stochastic transition instance t j ∈ T stoch , whereby H s = {h st |h st : R
-G : T → EXP is a guard function that assigns to each transition t ∈ T a guard expression of the Boolean type. m 0 : P → EXP is an initialisation function that assigns to each place p ∈ P an initialisation expression of a multiset type C(p) M S .
Here, R + 0 denotes the set of non-negative real numbers, • t j denotes the set of pre-places of a transition t j .
The formal semantics of GHPN C is defined by unfolding the coloured hybrid Petri nets into the equivalent low level one. Thus we firstly discuss how GHPN C can be unfolded into GHPN . Then, we extend the formal semantics which has been presented in [18] and apply it to GHPN C .
Unfolding Uncoloured Petri nets can be folded into coloured Petri nets, if partitions of the place and transition sets are given. Vice versa, coloured Petri nets with finite colour sets can be automatically unfolded into uncoloured Petri nets, which then allows the use of all of the existing powerful standard Petri net analysis techniques. The conversion between uncoloured and coloured Petri nets changes the style of representation, but does not change the actual net structure of the underlying biological reaction network.
In [32], we have presented an efficient unfolding method for coloured Petri nets with finite coloured sets, in which we provide two approaches to compute transition instances. For a transition, if the colour set of each variable in its guard has a finite integer domain, a constraint satisfaction approach is used to obtain all valid transition instances. Otherwise, a general unfolding algorithm is adopted, in which some optimisation techniques are used like partial binding-partial test and pattern matching [27,32].
Hybrid semantics The semantics of GHPN C is given in terms of the discrete and continuous transitions. To harmonise the mathematical notations, let T C = T − T D = T cont denote the set of continuous transitions. Before we proceed, we assume that the coloured expression which is assigned to an arc is evaluated to a numeric value. Likewise the initial marking.
Definition 2 (Enabling condition) Let N =< P, T, A, , C, F, V, G, m 0 > be a coloured generalised hybrid Petri net and m be the marking of N at time τ . A transition t j ∈ T is enabled in the marking m, denoted by m[t j , iff ∀p i ∈ • t j :
Definition 3 (Firing rule of discrete transitions) Let N =< P, T, A, , C, F, V, G, m 0 > be a coloured generalised hybrid Petri net, m a marking of N , and t j ∈ T D a transition enabled in the marking m, m[t j , at time τ . The transition t j can fire and reach a new marking m , denoted by m[t j m , at time τ + d j if it is still enabled at that new time, with:
where
is a delay which is associated to the discrete transition t j and d stoch (t j ) is the random firing delay with negative exponential probability density function calculated for each stochastic transition t j using its rate g(t j ).
According to the above enabling and firing definitions, discrete transitions follow a policy which is called an enabling memory policy [24].
Moreover, when multiple immediate transitions are concurrently enabled, conflicts are solved by computing the relative firing frequencies of each enabled immediate transition. That is if an immediate transition t j is enabled in the current marking m, then it fires with the probability given by (1).
where w is the weight assigned to immediate transitions.
The semantics of continuous transitions is analogue to the ones in continuous Petri nets with maximal firing speeds depending on time as introduced in [6] and tailored to the specific needs in systems biology. The transitions' current firing rates (instantaneous firing speeds) depend on the current marking of their pre-places (i.e., species concentrations). In what follows, the firing semantics of continuous transitions is formally given.
We introduce the following notation. Let v j (τ ) represent the current firing rate of a transition t j ∈ T C at time τ , m i (τ ) = m(p i ) denotes the current marking of a place p i at time τ , and f j (τ )=f (t j ) denotes the maximal firing rate of a transition t j at time τ , then:
Equation ( 2) implies that a continuous transition can fire with its maximal rate if it is enabled or its rate will be zero otherwise.
When a continuous transition is enabled, it fires as soon as possible and its effect on the connected places can be given by the following definition. T, A, , C, F, V, G, m 0 > be a coloured generalised hybrid Petri net, m a marking of N , t j ∈ T C a transition enabled in the marking m, m[t j , at time τ , and v j (τ ) denotes the current firing rate of the transition t j . The transition t j fires with:
Equations ( 3) and ( 4) are called outflow and inflow of a place p i , respectively, due to the firing of a transition t j [6].
Generation of the corresponding ODEs For a given transition t j ∈ T C , the functions read(u, p i ), inhibit(u, p i ) are defined as follows:
Then the ODE corresponding to each continuous place in GHPN C can be generated using ( 5)
The semantics of GHPN C which we have discussed in Section 4.2 can be produced via simulation. In simulating GHPN C models, the following two main issues need be considered: the partitioning of the net's transitions into stochastic and continuous ones, and the synchronisation between the stochastic and continuous regimes.
The partitioning of the net elements has an important influence on both the accuracy and the efficiency of GHPN C simulation. The more continuous transitions we use, the faster the simulation speed we get. However, this will have a high impact on the result accuracy. Thus, it is important to appropriately partition the model transitions into stochastic and continuous ones. This can be done using either static or dynamic partitioning. In static partitioning, the process of assignment of each transition to either the stochastic or continuous category is done off line (i.e., before the simulation starts). Transitions with high firing rate will be considered continuously while transitions with low firing rates are considered stochastically. Additionally, the number of tokens of the transition's pre-places plays a role in determining the transition type. An important merit of this approach is that there is no additional overhead to the simulation due to the partitioning process. However, transition rates can dramatically change during the simulation. This motivates the use of dynamic partitioning. Using dynamic partitioning, transitions can change their type during simulation from stochastic to continuous or vice versa. Such change is usually based on either the transition rate and/or the number of tokens in the transition's pre-paces. Nevertheless, dynamic partitioning will result in additional overhead due to the repeated check of fulfilment of the partitioning conditions. Therefore, if the gain due to the use of dynamic partitioning is less than the partitioning overhead, static partitioning should be used instead. A detailed discussion about this issue can be found in [18].
Similarly, the synchronisation between stochastic and continuous regimes is vital for the accuracy of the simulation result. At the end, stochastic and continuous transitions are not isolated from each other. They mutually affect each other. Thus, we need a mechanism to better perform the synchronisation. One choice is to use (6) to detect the occurrence of a stochastic event while the numerical integration is progressing [12],
where ξ is a random number exponentially distributed with a unit mean, and a s 0 (x) is the cumulative rate of all stochastic transitions. Once we agreed on a selection of a partitioning algorithm and a synchronisation approach we can execute the GHPN C model by repeating a set of steps. In [18], we have presented a hybrid simulation algorithm based on (6) using both static and dynamic partitioning. The main steps are outlined here while the details can be found in [18].
Starting from an initial marking, the simulation algorithm computes state changes over the time which are regarded as the current marking at each simulation step. Initially, the rates of stochastic and continuous transitions are calculated. Next, the accumulated rate function of stochastic transitions is calculated. When the model contains one or more continuous transitions, the simulation algorithm constructs an ODE for each place. ODE construction is done via (5). Afterwords, the simulator numerically integrates the resulting set of ODEs simultaneously with (6) until an event occurs. The event can be one of the following: the enabling of an immediate transition, the enabling of a deterministically delayed transition, a deterministically delayed transition has finished its delay, or a stochastic event has occurred (i.e., Equation ( 6) is satisfied).
In all cases the ODE integrator is interrupted and the appropriate action is taken (e.g., by firing a discrete transition). After that the ODE integrator is restarted with the new marking to account for the discontinuity which is due to the firing of a discrete transition.
When there is no continuous transition in the system being simulated, the simulation of GHPN C model is simplified to the simulation of stochastic Petri nets which can be carried out using, e.g., Gillespie's stochastic simulation algorithm [13].
All features of GHPN C which have been discussed in this paper are implemented in Snoopy [15] -a unifying Petri net tool that supports the construction and execution of different Petri net classes among them are: standard, extended, continuous, stochastic, and hybrid Petri nets. Those classes are available both on the uncoloured and coloured level. Snoopy is platform-independent and provided free of charge for academic purposes.
Snoopy also supports the automatic unfolding of GHPN C . In Snoopy, we can explicitly unfold a coloured stochastic/continuous/hybrid Petri net to its counterpart, a stochastic/continuous/hybrid Petri net, or we can do this implicitly during the execution of the Petri net model.
Furthermore, we have recently released a server-based version of the Snoopy simulator called S 4 (Snoopy Steering and Simulation Server) [19,20]. S 4 also supports the simulation of GHPN C . Besides, an GHPN C model submitted to S 4 can be steered (i.e., key simulation parameters can be changed on the fly), and collaboratively executed by different users. Snoopy and S 4 can be downloaded from http://www-dssz.informatik.tu-cottbus.de/snoopy.html.
In this section we present one case study as a typical application of Coloured Generalised Hybrid Petri Nets.
We will demonstrate the GHPN C using an example of a synthetic circuitthe repressilator, which is an engineered synthetic system encoded on a plasmid, and designed to exhibit oscillations [7]. The repressilator system is a regulatory cycle of three genes, denoted by, e.g., gene a, gene b and gene c, where each gene represses its successor, namely, gene a inhibits gene b, gene b inhibits gene c, and gene c inhibits gene a. This negative regulation is realised by the repressors, protein a, protein b and protein c, generated by the genes gene a, gene b and gene c, respectively.
The 1-bounded places as determined by P-invariant analysis and the related transitions as determined by T-invariant analysis are kept discrete. The unbounded places and related transitions are approximated by continuous places and transitions, respectively. That is, places gene i and blocked i, and transitions 1.
For example, if we assign the rates in Table 1 to the hybrid model and consider them as stochastic or deterministic rates, depending on the transition type, hybrid simulation yields plots as illustrated in Figure 3.
Figure 4 gives an GHPN C for the repressilator model in Figure 2. A colour set GeneSet is defined with three colours, a, b and c, to distinguish the three genes. Each place gets assigned this colour set GeneSet. A multiset expression, 1'all()=1'a++1'b++1'c, is assigned to the place protein. In Figure 4, we define a variable x of GeneSet, which is used in arc expressions. The predecessor operator "-" in the arc expression −x returns the predecessor of x in an ordered finite Fig. 3. Plot of one hybrid simulation run for the repressilator. For rate functions, see Table 1. This plot suggests that GHPN are able to capture the oscillation. Repeated runs look differently; thus stochasticity is captured as well. The coloured Petri net model in Figure 4 when unfolded yields the same uncoloured Petri net model in Figure 2. That is, the GHPN C model reduces the size of the original stochastic Petri net model to one third, which is a big advantage of coloured Petri nets.
In this paper we have introduced a new class of coloured Petri nets called Coloured Generalised Hybrid Petri Nets (GHPN C ). GHPN C are particularly tailored to systems biologists' needs to model and analyse multiscales models. Moreover, GHPN C provide the interplay between stochastic and continuous regimes on the coloured level which eventually provide a tool to control both the accuracy and the speed of running simulation.
In [21] we have proved that marking-dependent arc weights are necessary to model certain biological systems. However, this feature is still not implemented in Snoopy for GHPN C . Thus we plan to add it in order to widen the types of models that can be constructed via GHPN C . Furthermore, the current simulation of GHPN C is carried out by firstly unfold the coloured model into an uncoloured one, then it is simulated on that level. Although the unfolding process is fully automated, it may take a considerable amount of time, particularly for bigger models. Thus if an error is observed at the (beginning) of the simulation due to modelling mistakes, the entire unfolding process will completely be repeated. Therefore, simulating GHPN C directly on the coloured level [9] will increase the productivity and eventually save a lot of time.
Finally, the presented case study in this paper aims to demonstrate how GHPN C works. More sophisticated biological models can be constructed in the future using GHPN C . For example, we are working on using GHPN C to model diffusion-reaction systems described by partial differential equations, where we first discretise the space into a set of grid cells, and then we obtain a GHPN C model by representing each grid cell as a colour [28].