Modeling and Simulation with Hybrid Functional Petri Nets of the Role of Interleukin-6 in Human Early Haematopoiesis Sylvie Troncale, Fariza Tahi, David Campard, Jean-Pierre Vannier, and Janine Guespin Pacific Symposium on Biocomputing 11:427-438(2006)


MODELING AND SIMULATION WITH HYBRID FUNCTIONAL PETRI NETS OF THE ROLE OF INTERLEUKIN-6 IN HUMAN EARLY HAEMATOPOIESIS
SYLVIE TRONCALE AND FARIZA TAHI LaMI, CNRS-UMR 8042, Universit´ d'Evry Val-d'Essonne, Genopole, France e DAVID CAMPARD AND JEAN-PIERRE VANNIER Laboratoire M.E.R.C.I (EA2122), Universit´ de Rouen, France e JANINE GUESPIN Laboratoire de Microbiologie de Rouen, France
The regulation of human haematop oiesis is a complex biological system with numerous interdep endent processes. In vivo Haematopoietic Stem Cells (HSCs) selfrenew so as to maintain a constant po ol of these cells. It would b e very interesting to maintain these cells in vitro, in view of their therap eutical imp ortance. Unfortunately, there is currently no known pro cess to activate HSCs self-renewal in vitro. Since the difficulties related to in vitro exp eriments, mo deling and simulating this pro cess is indisp ensable. Moreover, the complexity of haematop oiesis makes it necessary to integrate various functionalities: b oth discrete and continuous models as well as consumption and production of resources. We thus fo cus on the use of Hybrid Functional Petri Nets, which offer a numb er of features and flexibility. We b egin by mo deling and simulating the role of a sp ecific cytokine, interleukin-6, in the regulation of early haematop oiesis. Results obtained in silico lead to the disapp earence of HSCs, which is in agreement with in vitro results.

1

Intro duction

Haematopoiesis is a complex phenomenon leading to the continuous production of all types of mature blood cells. This process is ensured by a population of haematopoietic stem cells (HSCs) which are able by a process called "selfrenewal" to maintain a constant pool in vivo. Nevertheless, the prolonged renewal of stem cells is not reproducible in vitro. We then concentrated our work on the modeling and the simulation of this biological process. For this purpose, we needed a formalism which integrates a number of features since the complexity of the biological system. Thus, the notions of production and consumption of resources as well as the notions of discrete and continuous time must necessarily be taken into account. Among the formalisms existing in the literature, we focused on Hybrid Functional Petri Nets (HFPN) 1,2,3,4 . These nets can accumulate the notions of interest


and they also have the particularity of allowing the definition of functions. We expressed a hypothesis concerning the role of a specific cytokine, interleukin-6, in the fate of HSCs 5 . The receptor of this interleukin (sIL-6R) is thought to be a candidate for a positive feedback loop. A consequence of this hypothetical feedback circuit would be bi-stability of the HSCs phenotype. Thus, depending on the concentration of sIL-6R, HSCs either commit to differentiation to form mature blood cells irreversibly or they self-renew. In this way, IL-6 activation can be considered as an epigenetic switch 6 , responsible for bi-stability (self-renewal and differentiation). This hypothesis was modeled and simulated in silico before performing costly in vitro experiments. To test this hypothesis, a pulse of sIL-6R was used to generate a transient modification of the environment, thus modifying the stable state of HSCs. This paper is organized as follows: the biological context is presented in section 2. Section 3 presents a brief state-of-the art on Hybrid Petri Net formalism. In section 4, we explain the different steps in the construction of our model. Finally, the last section presents the results. 2 Biological Context

The regulation of HSCs involves numerous growth factors. We focus on IL-6 and the molecules directly associated, since this signalling pathway is known to play a central role in stem cell biology 5 . Receptors involved in recognition of IL-6 are IL-6R and gp130. The assembly of the complete signalling receptor is sequential and hierarchically ordered. IL-6R binds IL-6 with a low affinity (Figure 1). IL-6R exists in two forms: a membrane-bound form, called mIL6R, and a freely-soluble-form, called sIL-6R. The soluble form of the receptor, sIL-6R, is secreted in an autocrine fashion by HSCs. IL-6/mIL-6R or IL-6/sIL-6R complexes are recruited by 2 gp130 subunit (catalytic receptor) and trigger activation of several intra-cellular pathways, particularly JAK (Janus kinase), which subsequently activates the STAT proteins (Signal transducer and Activators of Transcription) 5 . gp130/IL-6/IL6R complexes are endocytosed. The activation of JAK/STAT pathway leads to upregulation of gp130 8 and mIL-6R 7 in order to replenish HSC in IL-6 receptors. JAK/STAT pathway plays also a ma jor role in self-renewal of HSCs, and neo synthesis of receptors subunits to IL-6 should maintain sensitivity to this cytokine. We assumed that this up regulation of gp130 and mIL-6R is a significant mechanism supporting self-renewal. Production of sIL-6R is mediated not at a level of gene transcription but by a cleavage mechanism in-


Figure 1. Signal ling pathway of Interleukin-6. IL-6/IL-6R complexes and gp130/IL6/IL-6R complexes form what activates the JAK/STAT pathway. The STAT act as transcription factors and enable gp130 and mIL-6R synthesis as well as expression of SOCS proteins which inhibit the JAKs. A part of the mIL-6R is cleaved in sIL-6R. Finally, according to the level of sIL-6R, the p ermissive HSC can self-renew or differentiate.

volving ADAM-protease 7 . Finally, the JAK/STAT pathway activates SOCS proteins (Suppressor of Cytokine Signaling) 5 . These SOCS act as inhibitors of the JAK/STAT pathway. Gp130 stimulation is crucial for maintenance and proliferation of HSCs. We assume that the strength of gp130 activation determines cell fate, in accord with recent ligand/receptor signalling threshold theories 9 . Indeed, stem cell self-renewal is favored by maintenance of ligand/receptor complex level above a critical threshold; otherwise, below this threshold, the stem cells commit to differentiation. These data therefore suggest the presence of an autocrine loop but the functionality of the positive feedback circuit remains to be demonstrated. Before performing costly in vitro or in vivo experiments, it is reasonable to start with in silico experiments. Moreover, the study of such regulation is difficult to carry out in vitro. Indeed, experimentally HSCs are usually isolated using a specific surface marker such as CD133. However, the ensuing "CD133 population" is heterogeneous,


since it includes stem cells and progenitors. Progenitors committed in lineage differentiation represent the most important fraction: 95% of the population. We shall call this subpopulation C (committed)-cells. The remaining 5% are composed of quiescent HSCs residing in G0 for 75% and recruited in the cellcycle for 25% 10 . Quiescent HSCs will be referred to as Pq -cells and primitive cycling cells as Pp -cells. The Pp subpopulation is assumed to be permissive to signals from the micro-environment. This intermediate state permits us to link primitive- to committed HSCs and to integrate signals that determine fate of primitive HSCs (self-renewal or differentiation). Thus, it is difficult to specifically study the evolution of HSCs, since they are not experimentally identifiable. 3 Hybrid Functional Petri Nets

To model the regulation of early haematopoiesis, we needed to take into account the molecular degradation as well as cellular evolution. Therefore, the notion of consumption and production of resources was important. Moreover, regulation of haematopoiesis has two time scales. We thus needed discrete time to model cellular evolution, since cells are enumerated, as well as continuous time to model molecular interactions. Consequently, the use of hybrid modeling conciliating discrete and continuous time was necessary. Hybrid Petri Nets represent a good tool to model and simulate such a process. There are currently two sorts of Hybrid Petri Nets in use: those which offer a maximum of functionalities so as to facilitate modeling 11,3,4,12 and those whose goal is to develop methods of analysis and model verification 13,14,15 . To enable formal analysis, the second class does not offer numerous functionalities. Since the complexity of the regulation of haematopoiesis requires a maximum of functionalities, we used the first class of Hybrid nets. Hybrid Petri Nets (HPN) 11 , defined by David and Al la, are the first Petri nets which integrate both discrete and continuous time: - Discrete "places" represent entities which can be numbered (thanks to tokens). - Continuous "places" represent entities which can not be numbered. A real number (representing a concentration for example) is associated with each continuous place. - Discrete transitions fire after a delay of time dt. - Continuous transitions continuously fire at a rate v(t). - Normal arcs activate transitions by consuming resources. - Inhibitory arcs inhibit transitions and consume resources. - Test arcs activate transitions without consuming resources.


The HPN model was improved by Drath et al. to obtain Hybrid Dynamic Petri Nets (HDN) 12 . The HDN model allows definition of functions on continuous arcs and continuous transitions but the consumed quantity has to be equivalent to the produced quantity. Drath et al. added the ob ject oriented paradigm to HDN to obtain Hybrid Object Nets (HON) 12 . The main purpose was to encapsulate subnets within ob ject frames. Hybrid Functional Petri Nets (HFPN) have been recently developed 1,2,3,4 to integrate properties of HPN, HDN and HON. These nets allow the definition of functions on discrete and continuous arcs and transitions. The notions of functions on arcs and transitions seem very useful to model dynamic biological systems, since they enable establishment of links between places which are not directly related in the model. Another characteristic of this formalism is that the quantity of consumed resources can be different from the quantity of produced resources. Since the HFPN model integrates a maximum of functionalities, they were chosen for our pro ject. Throughout the article, we will use HFPN notations depicted in Figure 2.

Figure 2. Notations used in Hybrid Functional Petri Nets (HFPN). Discrete transitions fire after a delay of time Td, whereas continuous transitions continuously fire at a rate v(t). The HFPN mo del uses two kinds of places (discrete and continuous) and three kinds of arcs (normal arcs, inhibitive arcs and test arcs).

4

The use of HFPN to mo del the role of IL-6R in the regulation of early haematop oiesis

Haematopoiesis is a complex biological process comprising multiple and interdependent regulation networks. We therefore built the model in several steps. First, we built a sub-model for cellular dynamics and a sub-model for


molecular dynamics. Then, we connected the two sub-models by modeling the autocrine secretion of sIL-6R by stem cells as well as the kinetics of the signalling pathway. Experimental parameters were integrated in the final model. They came either from our experimental model or from literature. In the sequel, only few examples of these parameters will be done. 4.1 Cel lular sub-model

The cellular model represents the evolution of each cell lineage as a function of time. This sub-model is built in discrete time, since cells can be numbered. The sub-model is shown in Figure 3. It contains three discrete places representing the three biological entities: quiescent stem cells Pq, permissive stem cells Pp, and cells committed to differentiation C. Discrete transitions represent all the processes that allow a cell to change its state: - A Pq cell can enter the cell cycle to become a Pp cell, and inversely (respectively, transitions T1 and T2 ). A quiescent cell (Pq ) is consumed, whereas one permissive cell (Pp ) is produced, and inversely. - A Pp cell can self-renew with symmetric division when the two daughter cells remain primitive (T3). After mitosis assumed a one day duration (delay of 1 on the transition T3 ), a permissive cell is consumed to produce two primitive daughter cells. - A Pp cell can self-renew with asymmetric division when the two daughter cells are different; one stays primitive and the other differentiates (T4). - A Pp cell can differentiate without division (T5). - A Pp cell can divide (T6) to form two committed daughter cells.

Figure 3. Cel lular model. Discrete places (circles) represent the three types of cells: Pq are quiescent cells, Pp are p ermissive cells, and C are committed cells. Discrete transitions contain a time delay after which transitions fire. The number of consumed and pro duced cells is indicated on the arcs.


4.2

Molecular sub-model

The second sub-model is the molecular model representing interactions between cytokines. It represents the association and the dissociation of the complexes described above. Since we use growth factor concentrations, and since the formation and dissociation of complexes are continuous phenomena in cells, this sub-model was built in continuous time (Figure 4).

Figure 4. Molecular model. Continuous places mo del the molecules involved in the regulation of haematopoiesis by IL-6, and continuous transitions represent biological processes such as asso ciation of the complexes. The transitions T7 and T8 mo del the association of the IL-6/IL-6R complexes whatever the form of the receptor (sIL-6R or mIL-6R), whereas transitions T11 and T12 symb olize their dissociation. In the same way, transitions T9 and T10 mo del the formation of the gp130/IL-6/IL-6R complexes, whereas transition T13 represents their disso ciation.

Each cytokine as well as each complex of interest is modeled by a continuous place. We have places for the molecules IL-6, sIL-6R, mIL-6R, gp130 and places for all complexes: IL-6/IL-6R and gp130/IL-6/IL-6R with the soluble or membrane receptor. Continuous transitions used in the model represent biological processes such as formation and dissociation of all the complexes.


4.3

Autocrine secretion of the receptor of interleukin-6 and signal ling pathway

Once both sub-models were built (the cellular sub-model and the molecular sub-model), we modeled the influence of IL-6R on cells. We thus numbered the biological links between molecular interactions and the cellular evolution. The model is presented on Figure 5. As explained above, HSCs constitutively secrete the soluble receptor of IL-6, sIL-6R. In order to take this into account, we added a continuous transition (since Pp cells secrete a certain amount of molecules) where the entering arc is a test arc (since permissive cells are not consumed when secreting the receptors) (transition T14 on Figure 5). The in silico quantity of sIL-6R secreted by permissive cells was determined thanks to 7 and 16 , where the sIL-6R basal production by Pp cells was experimentally estimated to 43pg/mL/24h/150000 Pp -cells. The most important link between cellular evolution and molecular interaction is the kinetics of the signalling pathway previously described. When the complex gp130/IL-6/IL-6R is formed, it activates the signal for the activation of JAK/STAT pathway. After one hour, all the gp130 are endocytosed, which interrupts signal activation. We therefore used a discrete transition with a delay of one hour to represent this process (T17 on Figure 5). During this delay, gp130/IL-6/IL-6R complexes are not consumed, which is why we used a test arc. At the end of this delay, the virtual place Pv 1 (without biological meaning) is filled. Consequently, the T18 transition can fire triggering consumption of all the complexes. According to the quantity of formed complexes, either self-renewal is activated (Pv 2 and Pv 3 ) or cells commit to differentiation (Pv 4). For example, if the quantity of gp130/IL-6/IL-6R is sufficient, the virtual place Pv 2 is filled; the transition T3 can then fire and Pp -cells self-renew. The JAK/STAT pathway is activated as soon as the concentration of gp130/IL-6/IL-6R is above a certain threshold. Then, the SOCS proteins inhibit this pathway. This negative loop introduces a refractory period, chosen as 8 hours, in analogy with 17 . After this delay (modeled by the discrete transition T19 ), the gp130 and newly synthesized mIL-6R return to the membrane. So the cell becomes once again responsive to the formation of complexes. This receptor accessibility is symbolized by the T20 transition, which increases membrane receptor concentration. Finally, thanks to the JAK pathway, a part of the synthesized mIL-6R is cleaved into soluble receptor sIL-6R. The cleavage is represented by the T21 transition.


Figure 5. Global model. The cellular and molecular sub-models were joined thanks to biological links: the autocrine secretion of sIL-6R by p ermissive cells (T14) and the kinetics of signalling pathway (Pv 1-Pv 6, T17-T20). The hypothetical retroaction lo op is indicated in boldface arcs.

5 5.1

Simulation results and discussion Results obtained with the model

Once the appropriate model of the IL-6-mediated regulation of the early haematopoiesis was built, a set of simulations was carried out with exper-


imental values (detailed in 7 ). The simulations were done using the Cel l Il lustrator software 3,4 which implements the HFPN formalism. We concentrate our work on the study of the evolution of permissive cells (Pp ) as a function of time. The simulation of the model presented in Figure 6 leads to the disappearance of permissive cells in about ten days. This experiment was subsequently tested in vitro 18 , and the obtained results conform to the in silico results. Consequently, our model is adequate for further testing our hypothesis of an epigenetic switch.

Figure 6. Evolution of Pp cel ls as a function of time.

5.2

Test of an epigenetic modification

To test the functionality of the feedback circuit (in bold in Figure 5), it is necessary to verify whether a transient modification (using a pulse) of the sIL-6R (the candidate for the loop) is sufficient to activate the retrocontrol loop, and then to change the phenotype of permissive cells (self-renewal of stem cells). A pulse consists of adding a large quantity of a molecule, and after a delay, this molecule is completely removed. If haematopoietic stem cells are sub ject to an epigenetic modification, the feedback circuit would stay activated thanks to the inductor (sIL-6R), which exerts a positive regulation on its own synthesis. Two kinds of pulse were modeled: a pulse of sIL-6R and a pulse of HIL-6 (a fusion protein equivalent to the complex IL-6/IL-6R). Figure 7 presents the evolution of the number of Pp cells as a function of time, when a two hour pulse was applied. The symmetric self-renewal was activated, but only for one cycle after the pulse. Indeed, we observed a doubling of the totality of Pp cells. However, self-renewal was not activated longer and permissive cells became committed to differentiation.


Figure 7. Evolution of Pp cel ls as a function of time, when sIL-6R and HIL-6 are pulsed (confounded curves).

Simulation results revealed that there is no significant difference between the two kinds of pulses (sIL-6R and HIL-6). We thus conclude that the strong dissociation coefficient (K d = 10-9 ) of the complex IL-6/IL-6R do not influence the activation of the loop. Also, all permissive cells committed to differentiation despite the pulses. These results suggest that a pulse of sIL-6R (or HIL-6) with the experimental parameters is not sufficient to maintain stem cells in silico. Finally, the results of the simulations led us to suggest that activated cells produce sIL-6R in a quantity smaller than during the precedent cycle. Thus, each day a smaller fraction of Pp cells can self-renew (peaks observed every day). Consequently, a progressively larger part of Pp cells commits to differentiation. Note that these experiments were also performed in vitro 18 and the results confirm the results obtained in silico. 6 Conclusion and future work

Haematopoiesis is a complex process involving numerous interdependent circuits in its regulation. A formalism which integrates a maximum of functionalities is therefore needed. The Hybrid Functional Petri Nets present a high level of integration. They enable the modeling of production and consumption of resources at a time scale which can be discrete or continuous. They also allow the definition of functions on arcs and transitions. Thanks to the numerous functionalities available with this formalism, we successed in modeling the entire process of early haematopoiesis by interleukin-6. The simulation results predicted the disappearance of primitive HSC sub-


population whatever the tested pulse. This result suggests that this positive loop is not functional and that production of sIL-6R does not constitute an epigenetic modification determining the fate of HSCs, in the actual experimental conditions. The agreement of in silico and in vitro results provides informal validation of our model. Nevertheless, it would be interesting to validate a model in a more formal manner. At present, the HFPN model only enables simulations. We are therefore interested in developing a system of validation and verification for hybrid petri nets considering a maximum of functionalities. Finally, on a practical level, it would be useful in future work to vary those parameters whose values were only estimated (from the literature), i.e., those that could not be measured in our experimental model 7 . Also, we plan to add other factors, such as TGF- (Transforming Growth Factor) which maintains primitive cells in quiescent state. References 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. Nagasaki M et al, Applied Bioinformatics 2, 181-184 (2003) Doi A et al, Applied Bioinformatics 2(3), 185-188 (2003) Matsuno H et al , In Silico Biology 3, 389-404 (2003) Doi A et al, In Silico Biology 4, 271-291 (2004) Heinrich PC et al, Biochem J 374, 1-20 (2003) Thomas R and d'Ari R, CRC Press , (1990) Campard D et al, Stem Cel ls (submitted) , (2005) O'Brien CA and Manolagas SC, 272, 15003-15010 (1997) Viswanathan S et al, Stem Cel ls 20, 119-138 (2002) Cheshier SH et al, Proc Natl Acad Sci USA 96, 3120-3125 (1999) David R, Alla H and Bail J, Proc. 1st Int. ECC European Control Conference, Grenoble , 1472-1477 (1991) Drath R, 3rd International Conference on Automation of Mixed Processes, ADPM'98 , (1998) Wolter K and Zisowsky A, In Proceedings of the 4th International Computer Performance and Dependability Symposium IPDS'00 , (2000) Muller Ch and Rake H, Symposium on Mathematical Model ling , 457-460 ¨ (2000) J¨rns C and Litz L, In Proceeding of the 3rd European Control Conference o ECC95 , 2035-2040 (1995) Grafte-Faure S et al, Br J Haematol 105, 33-39 (1999) Gerhartz C et al, Eur. J. Biochem 223, 265-274 (1994) G¨tze K.S et al, Experimental Heamtology 29, 822-832 (2002) o