On determining firing delay time of transitions for Petri net based signaling pathways by introducing stochastic decision rules

---

Abstract

Parameter determination is important in modeling and simulating biological pathways including signaling pathways. Parameters are determined according to biological facts obtained from biological experiments and scientific publications. However, such reliable data describing detailed reactions are not reported in most cases. This prompted us to develop a general methodology of determining the parameters of a model in the case of that no information of the underlying biological facts is provided. In this study, we use the Petri net approach for modeling signaling pathways, and propose a method to determine firing delay times of transitions for Petri net models of signaling pathways by introducing stochastic decision rules. Petri net technology provides a powerful approach to modeling and simulating various concurrent systems, and recently have been widely accepted as a description method for biological pathways. Our method enables to determine the range of firing delay time which realizes smooth token flows in the Petri net model of a signaling pathway. The availability of this method has been confirmed by the results of an application to the interleukin-1 induced signaling pathway.

Keywords: Petri net, interleukin-1 (IL-1) signaling pathway, firing delay time, stochastic decision rule, conflict resolution

---

Introduction

A Petri net is a formal description for modeling concurrent systems [1]. Petri nets have recently been widely accepted as a description method for biological pathways such as gene regulation networks, metabolic pathways and signaling pathways by researchers in computer science as well as those in biochemistry [2]. Various types of Petri nets (e. g. stochastic Petri nets [3, 4], hybrid Petri nets [5, 6], colored Petri nets [7]) have been applied to study biological pathways in both quantitative and qualitative approaches because of potential advantages of Petri nets providing an intuitive graphical representation and capabilities for mathematical analysis.

Signaling pathways regulate elaborated cell communication mechanisms by controlling various alteration procedures of cell behavior such as cell growth, survival, proliferation, and apoptosis. By such cellular communication mechanisms, cell activities could be subtly governed and maintained in a good condition along with other biochemical interactions and processes.

Recently, Li et al. [8] proposed a qualitative modeling method of signaling pathways by paying attention to the molecular interactions and mechanisms using discrete Petri nets. Further, they proposed a timed Petri net based method of determining the firing delay time (or delay time for short) of transitions, i. e., the time each transition takes to fire [9]. Their method to determine delay times in a timed Petri net is based on the assumption that, for any reaction, a total amount of consumed substrates is equal to a total amount of products of the reaction. Although this assumption ensures the concentration equilibrium of the reaction, their method can only produce "exact delay time" for any transition in the Petri net model, i. e., no time range is allowed for any transition in determining delay time. This means that their method eliminates the variety of possible reaction speeds, which commonly exists in an ordinary signaling pathway. Furthermore, their method is designed based on the strategy that in the conflict situation at a certain place, "the same firing frequencies" should be assigned to the transitions going out from that place. This strategy should be improved, since it does not reflect real reactions in a cell.

This paper proposes a new method to determine delay times in a signaling pathway, which resolves above two problems, "exact delay time" and "the same firing frequency," while keeping smooth signal flows of a signaling pathway. A new concept "retention-free" is introduced to the Petri net. In a retention-free Petri net, the total token amount flowed out of a place per time unit should not be less than the one that flowed into the place. This new concept allows us to resolve the former problem, providing a delay time in some range to any of the transitions. On the other hand, the latter problem is resolved by introducing stochastic factors into transitions coming from a certain place where a conflict happens. These stochastic factors reflect the reaction rates of corresponding biological reactions, making it possible to realize more realistic signal flows of a signaling pathway.

The organization of this paper is as follows. After giving necessary definitions of timed Petri nets for the discussion of this paper, a Petri net modeling method of a signaling pathway is introduced with the example of the IL-1 signaling pathway according to the method proposed by Li et al. [8]. We then present formulas in two cases, conflict-free transitions and conflict transitions, which express the conditions for smooth token flows in a timed Petri net, i. e., smooth signal transductions in a signaling pathway. Based on these formulas, we propose an algorithm to determine the delay time of transitions in a timed Petri net model of a signaling pathway. The availability of this algorithm was confirmed by the results of applying this algorithm to the timed Petri model of the IL-1 signaling pathway.

---

Definitions of Petri nets

Petri net technology provides a powerful approach to modeling and simulating various concurrent systems [1], which has been widely employed as a description method for biologists and computational biologists owing to following advantages [10]:

1. "firm mathematical foundation" enabling formal and clear description of biological pathways as well as their structural analysis, and
2. "visual representation of networks" which provides intuitive understanding of biological pathways without any mathematical descriptions that are basically difficult for ordinary biologists.

We here briefly give the necessary definitions of Petri net and its extension (i. e., timed Petri net) used in the paper. For detailed definitions the reader is referred to [1].

Basic definitions of Petri net

A Petri net is comprised of three types of elements &ndash places, transitions and arcs &ndash whose symbols are illustrated in Fig. 1.

Definition 1: A Petri net is denoted as a 5-tuple PN = (T, P, E, α, β) that is a bipartite graph, where E = E⁺ ∪ E⁻,
T: a set of transitions {t₁, t₂, ..., t_|T|}
P: a set of places {p₁, p₂, ..., p_|T|}
E⁺: a set of arcs from transitions to places e=(t,p)
E⁻: a set of arcs from places to transitions e=(p,t)
α: α(e) is the weight of arc e=(p, t)
β: β(e) is the weight of arc e=(t, p)

Definition 2: Let PN be a Petri net.

1. •t (or t•) is a set of input (or output) places of t, and •p (or p•) is a set of input (or output) transitions of p.
2. A transition without input arc (Fig. 2) is called source transition. Let T_sour = {t₁^sour, ..., t_a^sour} (a ≥ 1) be the set of such source transitions. A source transition is always firable.
3. A transition without output arc (Fig. 2) is called sink transition. The set of such sink transitions is denoted by T_sink = {t₁^sink, ..., t_b^sink} (b ≥ 1).
4. A transition connected with two or more input arcs (Fig. 2) is called synchronous transition, and is defined by T_sync = {t₁^sync, ..., t_c^sync} (c ≥ 1).

A place can hold a positive integer that represents the number of tokens. An assignment of tokens in places expressed in form of a vector is called marking M, which varies during execution of a Petri net.

Firing rule of Petri net PN: A transition t is firable if each input place p_I of PN has at least α_e(e=(p_I, t)) tokens, where α_e denotes the weight of an arc e=(p_I, t). Firing of a transition t removes α_e tokens from each input place p_I of t and deposits β_e(e=(t, p₀)) tokens to each output place p_O of t. Fig. 3 shows the movement of tokens by the firing of the transition.

Conflict-free Petri net

A conflict (see Fig. 4) corresponds to the condition of a place that has at least two output transitions and has an insufficient number of tokens for firing. Deterministic resolutions are considered to decide the firings of the transitions. Usually, the concepts of priority-order and probability are employed when dealing with timed Petri nets, stochastic Petri nets and high-level Petri nets [11, 12]. A Petri net is called a conflict-free Petri net if no such conflict exists, which is a subclass of Petri nets.

Timed Petri net

a Petri net can be extended by assigning a delay time of firing to each transition for facilitating a system-level understanding by means of simulation. Such an extended Petri net is called a timed Petri net.

Definition 3: Let PN be a Petri net. A timed Petri net TPN is defined by TPN = (PN, D), where D is a set of firing delay times of each transition in T.

The firing rule of a timed Petri net TPN is defined as follows: (i) If the firing of a transition t_i is decided, tokens required for the firing are reserved. We call these tokens reserved tokens. (ii) When the delay time d_i of t_i passed, t_i fires to remove the reserved tokens from the input places of t_i and put non-reserved tokens into the output places of t_i. In a timed Petri net, firing times of a transition t_i per time unit is called firing frequency f_i. represents the maximum firing frequency of t_i. The delay time d_i of t_i is given as the reciprocal of .

Retention-free Petri net

Signaling pathways are composed of consecutive signaling events that are mediated by intracellular signaling proteins (usually enzymes) that relay the signal into the cell by activating the next enzyme from inactive state to active state on receipt of up-stream signal. With this feature, it can be considered that the signaling pathway will be in an abnormal state if the accumulation of the substances occurs per time unit. In Petri nets, such accumulation of a substance is represented by the token retention of the corresponding place. The token number of the place grows infinitely with the firing of its input transitions, and the lack of firing by the output transition. In this paper, such a Petri net in which signal flows can be steadily and smoothly propagated without the token retention at any place is called a retention-free Petri net.

Retention-free Petri nets are a subclass of timed Petri nets, where the total amount of in-flowing tokens () is not larger than the possible maximum number of out-flowing tokens () for each place p per time unit, which is represented by following inequality:

---

Modeling of signaling pathways

We give here a modeling method with a series of modeling rules for signaling pathways based on Petri net representation. As an example, the interleukin-1 (IL-1) signaling pathway is used to demonstrate our modeling method. Petri net based modeling of signaling pathways provides us with an intuitive understanding of the intrinsic structure and features of signaling pathways, and further enables computational experiments on the constructed Petri net model as being demonstrated in the following sections.

Modeling rules

The structural characteristics of signaling pathways can be naturally and explicitly expressed by Petri nets according to following rules.

1. Places denote static elements including chemical compounds, conditions, states, substances and cellular organelles participating in the biological pathways. Tokens indicate the presence of these elements. The number of tokens is given to represent the amount of chemical substances.
2. Transitions denote active elements including chemical reactions, events, actions, conversions and catalyzed reactions.
3. Directed arcs connecting the places and the transitions represent the relations between corresponding static elements and active elements. Arc weights α and β describe the quantities of substances required before and after a reaction, respectively. In case of modeling a chemical reaction, arc weights represent quantities given by stoichiometric equations of the reaction itself. Note that the weight of an arc is omitted if the weight is 1.

Signaling pathways are information cascades of enzyme reactions from transmembrane receptors to the nucleus DNA, which ultimately regulate intracellular responses such as programmed cellular proliferation, gene expression, differentiation, secretion and apoptosis. For signaling pathways, besides the catalytic reactions, the information among the molecular interactions such as complex formation, gathering action, translocation and channel switching need to be modeled according to different types of interactions as long as the biological facts are known.

Fig. 5 shows various molecular interactions of signaling pathways and their correspondence to the Petri net models. We give explanations about the molecular interactions that appear in the IL-1 signaling pathway, which is used below to demonstrate our method.

1. A binding reaction induces the formation of homo- or heterodimers and generates a complex compound. This block shows the ligand-receptor binding interaction and the corresponding Petri net model that indicates that the transition cannot fire in the absence of ligand although receptors exist. The number of input places of transitions is two or more while the output place number is one in the binding reaction. Obviously, we can also expand the conception of association to the model represented in block I (b), generally representing the simultaneous binding of substrates S₁, ..., S_n (n ≥ 1) forming a complex C in biological systems.
2. Phosphorylation is a reaction to add a phosphate (PO₃²⁻) group to a protein or a small molecule, and dephosphorylation that is the backward reaction removing phosphate groups from a compound by hydrolysis.
3. Autophosphorylation is a transphosphorylation reaction between receptor subunits frequently following the binding of a ligand to a receptor with intrinsic protein kinase activity.
4. Each down-regulated pathway transmits the signals to regulate different down-regulated reactions according to the modification positions on the ligand-receptor complex. The complex often has many chemical modifications, e. g., phosphorylation and acetylation. Few methods using Petri nets have been proposed to model this process by using one place (i. e., ligand-receptor complex) possessing more than one output transitions used to represent the multiple reactions of modification. These methods are easily understood, but raise an issue that, if the transition of such a place fires to remove the token(s) from a shared input place at one time, it will disable all the other transitions of this place simultaneously although the token(s) will return back to the same input place via a self-loop. To deal with it, we model each distinct active site (modification position) on the complex as an individual place C_i (1<i<n) as shown in block IV.
5. Adaptor protein (e. g., Grb2 and SHC1)-mediated association reactions are different from the binding reaction (see I). The main participator adaptor protein is an accessory protein to main proteins. These proteins lack the intrinsic enzymatic activities themselves but instead mediate specific protein-protein interactions driving the formation of protein complexes.
6. Chemical reactions are the most common reactions in signaling pathways, for which the conversion of substances to products is ordinarily modeled as connecting input places to output places, both belonging to the same transition.
7. Homodimerization is a polymerization reaction of two identical substances to shape a dimer similar to the association reaction. A substance is modeled as an input place connected with a 2-weighted arc. It is easy to expand the conception to model the formation of a multimer holding n-weight such as homotrimers or -tetramers.
8. Translocation refers to the movement of molecules, substances or ions across cell membranes or via the bloodstream in biology. Fig. 5 shows the nuclear translocation within a cell. A transition is modeled to indicate the movement action of substances before and after.
9. Intracellular signal pathways are largely carried out by second messenger molecules. Ca²⁺ acts as a second messenger molecule inside the cell. Usually the concentration of free Ca²⁺ within the cell is very low; it is stored inside of organelles, mostly the endoplasmic reticulum. In order to become active, Ca²⁺ has to be released from the organelles into the cytosol. Two transitions t_o and t_c are introduced to denote channel activity of "open" and "close", respectively. t_o is enabled when input place holds up token(s) after the association of organelles and substances, whereas t_c is enabled as long as some stop mechanisms shut off the channel.
10. This is the opposite of I. Dissociation process is a general process in which complexes and molecules separate or split into smaller molecules or ions. The number of input places of transitions is one while the output place number is two or more.
11. Since an enzyme itself plays a role as catalyst in biological pathways and there occurs no consumption in biochemical reactions, the reaction is modeled as a transition, where the enzyme is modeled as an enzyme place that has a self-loop with the same arc-weight. That is, once an enzyme place is occupied by a token, the token will return to the place again to keep the firable state, if the transition has fired.
12. A source transition represents an activity to provide substances that will take part in the reactions. A sink transition denotes small and natural degradation of a substance.
13. Internalization is a phenomenon to decrease the number of receptors on the surface of a cell membrane, due to be exposed to corresponding biological agents such as ligands for a long time. As a result of decreasing the number of receptors, responsiveness to the ligands is decreased. The internalization is modeled by an output transition connected with a place denoting the receptor via a normal arc.

IL-1 induced signaling pathway

In this paper, we use the example of the IL-1 induced signaling pathway (see Fig. 6) to demonstrate our proposed method. IL-1 is a proinflammatory cytokine, and plays an important role in regulating the mechanism of inflammation. the IL-1 signaling pathway is composed of the NF-κB (nuclear transcription factor-κB) pathway and the MAPK (mitogen activated protein kinase) pathway [13]. We construct the IL-1 signaling pathway by using our modeling method. Fig. 7 illustrates the Petri net model of the IL-1 signaling pathway.

---

Determining firing delay time of transitions

As shown in Fig. 7, the Petri net model without firing delay time of transitions is constructed. This model describes the structural information of connection relationships. The next task is to validate the model by means of simulation, i. e., to investigate whether the constructed model is consistent with biological facts obtained from biological experiments and scientific publications. Unfortunately, such reliable data describing detailed reactions are not reported in most cases. It thus leads us to develop a general methodology to determine the delay time of the transitions representing reaction rates.

Here we propose a new method of dealing with the delay time of transitions for a subclass of Petri nets, retention-free Petri nets, under two conditions. These are the conflict-free condition and the normal one including conflicts. We are to introduce a stochastic approach to determine the firings of transitions in conflict. Note that inhibitory arcs and cyclic structures (e. g., feedback loop, self-loop) are not taken into account in this contribution. That is, the Petri net dealt with here is an acyclic one without inhibitory arcs.

As defined above, in a times Petri net, the delay time d_i of t_i is the reciprocal of the maximum firing frequency , and it is obvious that the token amount flowed-in per time unit is no more than that flowed-out with . We thus decide d_i by calculating f_i.

Strategy for determining delay time in the case of conflict-free

In a conflict-free Petri net, each place may have one or more input arcs but at most one output arc. Without loss of generality, we are to consider two cases: (1) each transition has exact one input and one output arc and they construct a path from a source transition to a sink transition, as shown in Fig. 8; (2) There are l paths as in case (1) and these paths merge into a place p that is connected to a sink transition t, as shown in Fig. 9.

Delay times of transitions of case (1)
To keep retention-free behavior, each transition t_i on the path of Fig. 8 must satisfy , where is the token amount flowed into place p_i−1 per time unit and is the maximum firing frequency of t_i. It is obvious that is determined by firing frequency f_i−1 (not ) that is further determined by the firing frequencies of previous transitions. Since the source transition t₁ can fire times per time unit, that enables t₂ to fire times per time unit. In this way, Resultantly the following condition related to the maximum firing frequency of t_i can be obtained:

(1)

where is the maximum firing frequency of the source transition t₁, is the maximum firing frequency of t_i that is an output transition of place p_i−1, β_j (j=1, ..., i−1) and α_j (j=2, ..., i) are the weights of arcs connected from the input transitions and to the out transition as shown in Fig. 8.

Delay times of transitions of case (2)
For the transitions on each path of Fig. 8, the delay times can be determined according to inequality (1). And the token amount flowed into place p can be computed as discussed above and the maximum firing frequency of t can be determined according to the following inequality:

(2)

The left-hand of inequality (2) is designed to calculate total token amount flowed into p from its multiple input transitions per time unit.

Strategy for determining delay time in the case of conflict

Here, we propose a new firing rule of a timed Petri net by introducing a stochastic approach to determine firings for a series of transitions in conflict. Suppose a place p possesses output transitions, t₁, t₂, ... t_k then the firing rule is as follows.

New firing rule of timed Petri nets TPN*:

1. Each unreserved token deposited to input place p is assigned to be reserved by the transition t_i that satisfies the following mathematical expression:

   (3)

2. When the number of reserved tokens of t_i is not less than the required token number for the firing, the firing of t_i is decided.
3. After the delay time d_i of t_i passed, t_i fires to remove the reserved tokens from the input place of t_i and deposit unreserved tokens into the output places of t_i.

In the above expression (3), α_i is the arc weight of e(p, t_i); s_i is the firing probability of transition t_i, which represents the proportion of the firing frequency of each transition in the total firing frequency of the transitions in conflict. As shown in Fig. 10, s_i is assigned to corresponding transition t_i, which is given as a constant in advance according to the biological facts provided in the published literature. n_i is used to record the number of tokens that t_i has reserved so far; represents the total firing number of transition t_i from the beginning. Expression (3) is designed to reserve the token to such a transition t_i that has the largest difference between calculated firing probabilities and given firing probability s_i among all the transitions in conflict.

Delay times of transitions in conflict
We consider the determination of delay times for the transitions in conflict as shown in Fig. 11. The maximum firing frequency must satisfy the following inequality according to the above firing rule:

(4)

where α_j and β_i are the weights of and , respectively; s_j is the firing probability of ; is the firing frequency and is the maximum firing frequency of and , respectively, as shown in Fig. 11. The left-hand of expression (4) derives the token amount reserved to the transition from the total token amount coming from each input transition connected to the place p_i, while the right-hand represents the possible token amount from p_i to its output transition per time unit. represents the ratio of the token amount deposited to to the total token amount from p_i to each transition per time unit.

A synchronous transition t_i^sync connected from two or more input places becomes firable if all the input places have a sufficient amount of tokens satisfying the firing condition. Due to the feature of retention-free Petri net, it is necessary to synchronize the timing that each input place of t_i^sync becomes firable. It leads us to consider the condition to determine the delay time for the synchronous transitions in order to avoid token-retention at multiple input places merged to the same synchronous transition. The condition is, if there is a synchronous transition t_o^sync as shown Fig. 12, the maximum firing frequency satisfies the following equation:

(5)

Equation (5) is an expression formulated by the firing frequency of each input place of the synchronous transition obtained by using expressions (1), (2) and (4).

Algorithm of automatically determining delay time of each transition

We develop an algorithm to automatically determine delay time of transitions according to expressions (1), (2) and (4). First we give the explanation of the notions used in our algorithm.

LT is a list, , containing synchronous and sink transitions in such order that if i<j then there is no directed path from t_j^sync to t_i^sync. T_sync and T_sour are sets of synchronous and source transitions, respectively. EQ1(t^sync) is a function that outputs conditional expression of a synchronous transition t^sync based on equation (5). EQ2(p) is to produce expression (4) for a conflict place p and EQ3(p) is to produce expressions (1) and (2) for a conflict-free place p. push(stack,p) is a stack operation pushing a place p into the stack stack. pop(stack) is a stack operation pulling out from stack. T_e is a set of transitions whose delay times have been derived from EQ2(p) and EQ3(p).

In the above algorithm, the main routine Eq-Produce(PN) is designed to find all the synchronous transitions for delay time determination. In Step 1º, synchronous transitions are sequentially saved to the list LT from the source transitions to the sink transitions. The sets of synchronous and sink transitions are initialized in 2º and 3º, respectively. In 5º, the sub-routine dfs-push(t, stack, PN) is executed, in which first argument t is a transition pulled from the head of list LT in 4º. In 6º, the function EQ1(t^sync) outputs conditional expression regarding t^sync by using equation (5). The procedures in 4º - 8º are executed until LT is empty. That is, if all the synchronous transitions stored in LT are searched, the algorithm will stop.

The sub-routine dfs-push(t, stack, PN) is a function recursively carrying out an operation of Depth First Search from selected transition t (given as the first argument in the main function) towards the source transitions. In 2º, a judgment procedure if the input place p of t is a conflict place or not, is executed. When p is a conflict state, it will be labeled and pushed to stack. Steps 4º - 7º are the procedures for the input transitions t of p. If the condition t ∈ T_sour ∪ T_sync ∪ T_e holds, dfs-push(t, stack, PN) will stop and the other sub-routine dfs-pop(stack, PN) will be invoked; otherwise, dfs-push(t, stack, PN) will be invoked recursively.

The sub-routine dfs-push(t, stack, PN) is a function used to produce conditional expression with the use of places saved in stack based on expressions (1), (2) and (4). In 3º, if the input transitions of p are Marked, p will be pulled out from stack. In the case that p is a conflict place, EQ2(p) is used to calculate conditional expression (in 4º); in the case that p is a normal place, EQ3(p) is used (in 5º).

Application of proposed algorithm to IL-1 signaling pathway

As a case study, we applied our algorithm to a IL-1 signaling pathway model. Here, we show the results by applying our algorithm to automatically produce delay times of all the transitions to the IL-1 signaling pathway model shown in Fig. 7. In the IL-1 signaling pathway, it can be found that there exist several self-loops. We thus break down the IL-1 signaling pathway at the places of self-loops, and derive all the delay times of the transitions by applying the algorithm <<Firing Time Determination>>. Calculated delay times of the transitions are shown in Tab. 1, which are given in the form of conditional expression of firing frequency. The simulation can be executed without token-retention by using the delay time decided obeying the conditional expressions listed in Tab. 1. If the simulation is conducted in sufficient time period, it can be observed that firings of the transitions follow the given stochastic values.

---

Conclusions

We have proposed a method of determining the delay time of transitions satisfying derived conditional expressions for a class of Petri nets, which is acyclic and of no inhibitory arcs. To resolve nondeterministic firings, we have introduced a stochastic approach to determine firings of transitions in conflict.

In this contribution, we have first presented basic definitions of Petri net and introduced a Petri net based modeling method for signaling pathways by taking notice on the molecular interactions and mechanisms. Then we have presented conditional expressions for two cases, conflict-free transitions and conflict transitions, which express the conditions for smooth token flows in a timed Petri net, i. e., smooth signal transductions in a signaling pathway. Based on these expressions, we have proposed an algorithm to determine the delay times of transitions in a timed Petri net model of a signaling pathway. The availability of the algorithm has been confirmed by the results of application to the IL-1 signaling pathway model.

Using our method, the range of delay times of transitions can be decided by introducing the concept of probability to the conflict transitions. Our method makes it possible to automatically determine the delay times of all the transitions representing biological reaction rates in signaling pathways according to the delay times of reliably determined transitions. In the meanwhile, when several delay times of the transitions have been given in advance based on confirmed biological facts, the delay time of other transitions can also be determined mechanically.

The Petri net dealt with in this paper is constrained to an acyclic one without inhibitor arcs. However, it is obvious that there exist cyclic structures in signaling pathways such as feed-back loops and self-loops. Therefore, as the future work, we will aim to develop our current method further to deal with the Petri net model including various kinds of net structures.

---

Acknowledgements

This work was partially supported by Grant-in-Aid for Scientific Research on Priority Areas "Systems Genomics" (17017008) from Ministry of Education, Culture, Sports, Science and Technology of Japan.