BACKGROUND
the immune system  is a highly distributed system composed by cells and molecules whose scope is to protect living organisms against diseases.
the immune system has various defense mechanisms whose complexity is usually correlated to the host organism.
in mammals, among many defense mechanisms such as innate immune response arm, the adaptive is response, driven by lymphocytes, represents a fundamental component.
the goal is of the adaptive is response is to give the best reaction against most dangerous viruses and bacteria.
specificity, memory and discrimination between self and non-self represent the most important peculiarities of the adaptive is response.
in particular, the specificity feature is the capability, shared only by lymphocytes and antibodies, to recognize one specific epitope  of foreign pathogens  <cit> .
specific immune system lymphocytes can then recognize billions of different antigens, but for each determinant, a specific lymphocyte will be stimulated.
there are as many stimulated lymphocytes as determinants forming the antigen.
indeed, the memory feature represents the capability of the immune system to remember previously encountered antigens.
when an antigen appears for the first time, the immune system produces a first  response against it.
this stimulates the production of long-lived memory lymphocytes that are specific for each of the antigen epitopes.
when the same antigen comes again into contact with the immune system, the memory lymphocytes are stimulated in order to produce a faster and more effective immune system response against it .
building a model  of the immune system may lead to multiple benefits, since researchers can use immune system simulation in vaccine research, or to test biological hypotheses about diseases and infection dynamics.
many mathematical and computational techniques have been applied to model the immune system response; most of them bear on two main approaches: the top-down approach and the bottom-up approach.
if the former mainly focuses on the behavior of a system as a whole , the latter analyzes the system dynamics at a lower level  and the emergence of global properties is obtained as the sum of local behaviors.
models based on the top-down approach usually make use of differential equations like ordinary differential equations  .
these models are all population-based, and the spatiality and topology which both depend on individual interactions are, in general, ignored.
for very simple models it is possible to use some powerful mathematical techniques to extract some important system properties like analytical solutions, steady states and asymptotic behavior.
however, as the biological description introduced into the model grows, they become very complex and, in some cases, almost intractable.
indeed, the application of partial differential equations , stochastic differential equations  and delayed differential equations  can be introduced to obtain some spatial description, stochastic effects and temporal delays.
however, this usually entitles instability problems, a higher computational effort and numerical simulations only.
cellular automata and agent based models represent the most famous bottom-up approaches.
since entities are followed individually, they allow a much more rich description of the biological background, and individual cells properties like cell receptors, cell internal states and individual life-times and behaviors are not a problem.
furthermore, spatial description, stochastic effects and delays are intrinsic characteristics of these kind of models.
however, such techniques lack of tools and mathematical techniques for analytical studies and they usually require huge computational resources for large scale simulations.
petri nets  represent a graphical and mathematical modeling technique for the description of distributed systems presented in  <dig> by carl adam petri  <cit> .
during these last decades, petri net modeling formalism has been extended in order to manage, among others, stochastic effects, continuous quantities and scheduled events.
among pn extensions there is the colored petri nets  formalism  <cit> .
petri nets and extensions have been successfully used in bioinformatics for the modeling of intra-cellular biochemical and metabolic pathways .
park et al.
showed how it is possible to apply a particular extension of pn  for the modeling of the immune system response  <cit> .
the approach had several limitations; among them, it did not allow to discern between cells featuring different internal states, receptors, or a different position inside the simulation space.
consequently, such an approach did not bring any real additional advantage over up-down traditional techniques such as ode systems, and neither in respect to abms.
here we propose a novel methodological approach that uses “high level” petri nets, and in particular colored petri nets to model cell specific properties such as cell internal states and receptors.
to this end we present a proof-of-concept model of the humoral part of the is capable of reproducing the most important features of the adaptive response, and in particular the memory and specificity features.
the approach holds somewhere in the middle of top-down and bottom-up approaches, since it allows to describe cell specific features and behaviors, but with the possibility of using many standard pn analysis tools that permit to study system properties and features from a qualitative point of view.
all of this using a convenient graphical language that is already known to life scientists that, as already stated, are actually able to use pn for the modeling of signaling pathways.
the paper starts with an introduction to petri nets and the colored extension used in this paper.
details of the modeling approach are provided in methods, followed by results and discussion where we present the application of the approach to a real scenario, ending with our conclusions and final remarks.
methods
introduction to classical “low-level” petri nets
pns is a set of particular graphical symbols called: places, transitions, arcs and tokens.
the places are generally drawn as ellipses or circles and represent states and/or conditions of the system modeled; each place may be featured by a name , which has not a formal meaning but increases understanding of a petri net.
boxes or rectangular bars represent transitions which model events that may occur in the system.
places usually contain a discrete number of elements called tokens, drawn by black dots.
arcs may connect places to transitions and vice-versa, but they cannot connect transitions to transitions or places to places.
an arc connects a single  place to a single transition; furthermore, an arc may connect a transition to a single  place.
a transition may have several input and output places.
the tokens in the input places are fundamental for the firing  of the transition; the enabling of a transition models the occurrence of a specific event.
firing of a transition may require the presence of a different number of tokens for each input place.
generally, the required number of tokens for each input place is given by the multiplicity of each arc connecting input places with the transition.
when a transition fires , it consumes the required number of tokens from the input places, and creates tokens in its output places.
the number of tokens put in each output place depends again by the multiplicity associated to each arc.
a firing rule is atomic, in the sense that it represents a single non-interruptible action.
unless an execution policy is defined, when multiple transitions can be enabled at the same time anyone of them may fire .
introduction to “high-level” colored petri nets
a colored petri net  extend the definition of petri nets as it is a set of places, transitions, arcs and color sets.
each place has an associated type determining the kind of data that the place may contain .
like pns, during the execution of a cpn  each place will contain a varying number of tokens.
each of these tokens carries a data value that belongs to the type associated with the place.
for historical reasons, token values are referred as token colors and data types are referred as color sets.
this is a metaphoric picture where the tokens of “high-level” cpn are considered to be distinguishable from each other and hence “colored” – in contrast to ordinary “low-level” petri nets which have “black” indistinguishable tokens.
the types of a cpn can be arbitrarily complex e.g., a record where one field is a real, another a text string and a third a list of integers.
hence, it is much more adequate to imagine a continuum of colors  instead of a few discrete color values .
in cpns, there are different types of expressions.
an expression is mainly made up by variables.
it is not only associated with a particular color set but also written in terms of a predefined syntax.
in the following, we denote by exp a set of expressions that comply with a predefined syntax.
an expression may be associated to a transition; in this case, it is called guard expression.
an arc may have an expression associated, too; it is called arc expressions.
on the basis of what said until now, the formal definition of colored petri nets will be given in the following.
it uses the concept of multiset, which can be defined as a set in which there can be several occurrences for the same element.
the number of occurrences of an element is called coefficient or multiplicity.
for example, an infinite number of multisets exist which contain elements a and b, varying only by multiplicity, such as: {a, b}, {a, a, b}, {a, a, a, b, b, b}.
in the following, given a generic set e, the collection of all the multisets over e is denoted by e
ms.
definition: a colored petri net is a tuple n =  where:s is a finite non-empty set of types, called color sets.
p is a finite non-empty set of places.
t is a finite non-empty set of transitions.
a is a finite set of direct arcs
c is a color function that assigns to each place p∈p a color set c∈s.
c: p → s
g is the set of g functions that assigns to each transition t∈t a guard expression of the boolean type.
g: t → exp
f is the set of f function that assigns to each arc a∈a an arc expression of a multiset type cms, where p is the place connected to the arc a. f: a → exp
i is an initialization function that assigns to each place p∈p an expression of a multiset type cms.
i: p → exp
according to the tokens distribution in the set p∈p, the cpn may assume different states; each state is called marking and is given by a particular token distribution in each of the places p∈p.
the initial marking is produced by the application of the initialization function i. in the following the current token distribution in a particular place will be indicated by m, which is a function of a multiset type cms.
given a particular a transition t ∈ t, let us consider the set of variables present in the guard of t and in the arc expressions of arcs connected to t; this set will be denoted var.
before the guard and arc expressions relevant to the transition t ∈ t may be evaluated, the variables must get assigned values; this process is called binding.
a binding b of a transition t ∈ t is a function that maps each variable v ∈ var onto a value b that is of the same type as the variable.
the set of all bindings for a transition t is denoted b. a binding element is a pair , with t ∈ t and b ∈ b.
the behavior of the cpn is based on a firing rule consisting of a precondition and the effect of the occurrence  of a single transition.
whether or not a transition can fire depends on the marking of its input places and the arc expression on the input arcs.
a transition is enabled and is allowed to fire, if all the input places are sufficiently marked; this occurs when the binding of the variables that appear in the arc expressions of the arcs from inputs places to transition, evaluates to a multiset of token colors ms) that is present on the corresponding input place.
in this case, the guard of the transition should evaluate to true for the giving binding.
if a transition has no input places, it is always enabled.
in a more formal way, it is possible to say that the guard of the transition g ∈ g evaluates to true if m ≥ f, ∀p ∈ •t, where f ∈ f and •t is the set of input places for transition t. due to the definition of set f, f denotes the value ms) of the arc expression relevant to the arc going from p to t, for the giving binding b.
if the guard of the transition evaluates to true , the firing of a transition occurs.
in this case, a multiset of colored tokens are taken from each input place and added to output places in accordance with the values given by the arc expressions relevant to the arcs connecting the transition to output places.
this means that the firing process lead to a change in the current marking of the cpn.
in particular, the current marking m becomes the new marking m’, defined by m’ = m + f’-f, where f’ ∈ f is evaluated for ∀p ∈ t•  and denotes the values of the arc expressions relevant to the arcs going from p to t, for the giving binding b.
a binding element may occur concurrently to other binding elements – iff there are so many tokens that each binding element can get its “own share”.
the colored generalized stochastic petri net framework extension
different cpn modeling paradigms have been defined during years: colored qualitative petri nets , colored stochastic petri nets , and colored continuous petri nets   <cit> .
qpn
c does not uses tokens but it is only used for qualitative analysis of systems; spn
c takes into account discrete tokens and stochastic transitions, and cpn
c uses continuous quantities rather than discrete tokens and continuous transition rates instead of discrete transitions.
qpn
c is an abstraction of spn
c and cpn
c, while spn
c can reproduce by approximation cpn
c and vice-versa.
from our point of view, the most interesting approach is represented by spn
c. in brief, this colored pn extension uses a discrete number of tokens on its places like in classical pns .
moreover, stochastic events and transitions that are fundamental for simulating, for example, different individuals, are flawlessly reproduced through the use of stochastic transitions which fire after a probabilistic delay determined by a random variable defined according an exponential probability distribution.
for our purposes it is useful to consider a further extension of the spn
c framework, namely colored generalized stochastic petri nets, by considering, besides of stochastic transitions, also some kinds of deterministic transitions that we can identify as immediate, delayed, and scheduled transitions.
all transitions become enabled if all the preplaces are sufficiently marked.
when a stochastic transition  is enabled, a given time must be wait before the firing occurs.
this waiting time that determines the firing delay of the transition is given by a random variable x
t that is distributed exponentially with the following probability density function: fxtτ=λtm⋅e−λtm⋅τ,τ≥ <dig>
where λt can be defined as an arbitrary mathematical function which depends on the marking m of the preplaces at time t .
it should be noted that, even if there is a stochastic time delay before the firing of the transition, the firing itself does not consume any time.
deterministic  transitions  have a deterministic firing delay specified by an integer number.
the delay count starts just after the transition becomes enabled.
however, it must be said that during this waiting time it can happen that the transition loses its enabled state .
immediate transitions  can be seen as a special case of deterministic  transitions with a delay time set to  <dig>  in case of conflict between an immediate transition and any other kind of transition, the former will get firing priority.
also scheduled transitions  can be seen as a special case of deterministic transitions.
the firing is deterministic and it occurs at a previously defined absolute time of the simulation, obviously only if the transition is enabled at that time.
interested readers can find further details about spn
c  <cit> .
advantages of using colored petri nets
petri nets represent a graphical modeling tool that allows to describe in a simple and clear, but yet formally correct and powerful manner, any kind of process.
the biggest problem of classical low-level petri nets is given by the fact that they usually not scale.
as the  process we want to describe grows, low-level petri nets tend to grow quickly, and the designing, drawing, understanding and managing of the net becomes even more and more difficult, thus increasing the developing time and the risk of introducing modeling errors.
the introduction of “colors” allows to discern among different tokens and to follow the individually, in a similar way as it is done by ca and abm approaches.
this represents the most important feature of colored petri nets, and their greatest advantage in respect to low-level pns, but also in respect to all the differential equation based models.
just to give to the reader an idea of the practical advantage of cpns, the proof-of-concept model that we developed according to the approach described in the next section required only  <dig> places,  <dig> transitions and  <dig> arcs.
if a standard  pn approach would have been used, it would have required approximately  <dig> different places,  <dig>  transitions and  <dig>  arcs.
this means that the development of the model itself would have been almost unfeasible.
a similar scenario holds, for example, if a standard ode approach is used .
the biggest advantage of cpns  in respect to classical bottom up approaches is given by the availability of a large number of formal analysis methods by which the important properties of the net can be extrapolated, proved and analyzed.
liveness, boundedness, place invariants, circuits, t-invariants and so on  <cit> , allow to study the net properties from a qualitative point of view, making the pn approach more similar to a differential equation based approach in this sense.
furthermore, over classical simulation, it is possible to have interactive simulations with results drawn directly and in real time on the cpn diagram, and it is possible to act during these simulations by, for example, introducing and modifying new and existing tokens.
this kind of interactive simulation allows to debug the model behavior while it is developed, just like in a complete programming language environment.
indeed, novel simulation techniques that should speedup and potentially allow large scale simulations as, for example, simulations that use the “symbolic marking” of the net or “fluidification” techniques have now been presented  <cit>  and their implementation in present pn software is in due course.
finally, cpns  have a graphical representation.
this graphical representation is very simple, intuitive, and very appealing, even for people who are not very familiar with the formal details of cpn.
with only three kinds of graphical items  it is possible to mostly describe any kind of process and to obtain, for free, a very clear and qualitative conceptual model of process itself.
arc and transitions guards can include very simple first-order logic rules, as well as more complex functions that can be defined using look-up tables.
thus, nor a strong knowledge of mathematics , neither a good knowledge of any programming language  is required, and this may represent a consistent leap ahead in respect to classical approaches, since it could allow to life-scientists to directly use cpns with little or no effort at all.
how to use cpn to model the immune system function: a novel methodological approach
to develop the adaptive is response model, we used colored generalized stochastic petri nets, in which, as already described, the most common transition rules are defined as stochastic processes.
we used the snoopy software  <cit>  for developing a proof-of-concept model of the humoral response to an extracellular bacterium.
the approach we present can be reassumed as follows:use places to identify cellular or molecular types, or subpopulations types
use transitions to model events, such as cellular and molecular state changes, entity interactions, duplication, differentiation, aging and death.
use tokens to represent entity instances, like cells and molecules .
use colors to define specific properties of cellular and molecular entities.
in particular, we defined two distinct colorset types  that include the following colours: state, position, lifetime and receptor.
the state colorset is used to define internal cell states that determine cells behavior and it is actively used only for cellular types.
in the present implementation its possible values are: “resting; active; presenting; memory; plasma”.
the life colorset  is used to define life-time of cellular and molecular entities entities.
the receptor colorset can be used for many functions and its meaning depends on the cellular/molecular entity and transition rule we are dealing with.
it can represent:the cellular receptor that is used for the recognition of specific targets by b or t cells;
the presented peptide sequence that is presented by antigen presenting cells  like mp;
the epitope sequence presented bacteria or antigens in general;
the igg antibody specific binding structure.
the position colorset can be used to model position of cellular entities inside the space or to represent different compartments, i.e.
viral infection place and lymph node.
by using this approach and colored petri nets it will be possible to describe in a very general way, with just one transition, the behavior of all the entities that behave similarly , forgetting about the properties  that are not relevant for the behavior.
in other words it will be possible to say, for example, with just only one transition: “all the naive b cells that encounter in the same position  and recognize with sufficient affinity  a specific pathogen will be activated, no matter what is the valule of their real exact position, their value of their real exact receptor or peptide sequence, and how much remaining life they have”.
with this approach we will gain a more grained description that introduces, for example, cell receptors and allows to follow entities individually, but at the same time we will obtain a very easy but yet powerful and flexible way to define the general rules that drive the dynamics of the model.
of course, the definition of the used colors can be adapted to the problem we are dealing with i.e., the position colorset can be removed if we are not interested in modeling space and/or different compartments.
on the contrary, we can consider two or tree different position colors  to describe two or three dimensional environments.
in our proof-of-concept model we have taken into account  <dig> different places that are used to represent bacteria, macrophages , t helper cells , b cells , plasma cells , immunoglobulins class m  and immunoglobulins class g .
as already stated, we defined two distinct product colorsets  types for all tokens, one type that is mainly used for cellular entities and one type that is used for molecular entities.
this first product colorset is composed by four colors: state, position, lifetime and receptor, and it is used by th, b, pb, mp; the second one is composed by three colors: position, lifetime, and receptor and it is used for igm, igg, and bacteria.
we note here that, even if bacteria are cells and not molecules, we used the second colorset definition for such entities since we were not interested in describing bacteria internal states.
it is worth nothing that in the present model implementation all cells and molecules have distinct product colorsets  in order to improve readability.
however, it must be noted all the cellular entities have the same product colorset structure, and this also holds among molecular types.
therefore, it would be possible to define an unique product colorset for all cells and another one for all molecules, since the distinction among different cell types  or different molecular entities  depends on the context.
in fact, since places represent entity types, two identical tokens  may represent two different kinds of cells if they hold in different places.
the transition rules that are based on well-known biological rules define the correct state changes of the system and do not allow the erroneous migration of tokens towards places that represent different entity types.
it will be assumed that for any arc expression we can either have a variable with of a given product colorset type that will bind to any token with the same colorset, or a more complex expression, always based on the colorset definition that will select only tokens that will satisfy some given particular requirements.
for example the arc expression 2’, with p representing a variable whose type is position and l a variable whose type is life, when placed on an input arc  will only accept two tokens with subset colors state=“active” and receptor=“4”.
p and l will bind to any position and life colors of the incoming tokens.
in other words, any token representing a cell in the active state and receptor  <dig> will be accepted by the arc expression, no matter what life and position it will have.
the same arc expression when defined on an output arc  will produce two tokens with colors state=“active”, receptor=“4”, and position and life colors defined according to the p and l values.
if p and l binding already occurred, for example on the input arcs of the transition, their value will be already defined and assigned on the tokens that will be produced; if no binding of p and l occurred on the input arc expressions, their value will be chosen at random inside the set of the allowed position and life colorset values.
this example takes into account tokens which represent cellular entities, however, if we remove the state color set from the product color set definition, it is also valid for molecular entities.
as already stated, the evolution of the system is driven by transition rules that are defined upon well-known biological rules.
generally, transitions can be used to delineate all the biological scenarios that are needed to describe the evolution of the system.
among such scenarios we recall the appearing of new entities and death of old entities, cellular and molecular behavior changes that include changes in the internal cell state, aging , processing and presentation of antigenic sequences and movement.
such changes can be the result of transitions that involve one or more entities.
of course it is not possible to describe in detail all the possible biological scenarios that may arise, so we sketch out more explicitly the methodology by showing the most relevant examples.
introduction and death of new entities
the introduction and of new entities can be modeled as presented in fig.
<dig>  usually, new entities can be introduced with stochastic transitions with an initial pre-defined color.
for example, as shown in fig.
<dig>   <dig> th cells are introduced as “resting” cells, with an initial life-time set to  <dig>  both the receptor and the position colors are not assigned , and this translates into the generation of new cells with random position and receptor, but with a predefined life counter and initial internal state.
as it also possible to note from fig.
<dig>  such kind of transitions does not need any incoming token from a given place to be enabled.fig.
<dig> example of introduction of new entities into the system.
<dig> new th cells are introduced at random times into the system as resting cells , with a random position , life =  <dig>  and a random receptor
for some scenarios it is requested to model the introduction of new entities at a given time, for example to reproduce the injection of a given treatment or compound.
if such situation holds, it is possible to substitute the simple  transition with a scheduled transition that allows the definition of an initial time, period, and final time of activation of the transition.
death of entities will occur for all cells and/or molecules that have ended their life cycle .
this can be achieved by using transitions with transition guards on the value of the component element that refers to the life counter, as presented in the example shown in fig.
<dig>  here, any token can be potentially selected, in fact any token coming from the th place can bind to the variable th, but the transition will fire only if the guard condition  on the transition is satisfied, that is if the third component of th, representing the life subset color, is  <dig>  if such holds, the input token will be consumed by the transition, and no output tokens will be produced .fig.
<dig> example of the use of a stochastic transition to model the dead of entities that have completed their life cycle.
th cells that have completed their life cycle will be removed from the system.
the transition is enabled if and only if the guard condition on the life subset color holds, thus only if the binding on the th variable  is done on a token whose third component element  of the product colorset th is equal to zero
internal state changes
with the preposition “internal state changes”, we define all the events that modify the binding of the components of the product colorsets we have previously defined .
such changes may hold from an immunological point of view for multiple reasons.
a first example of this scenario can be represented by entities aging.
in this case the, life component of entities needs to be decreased to reflect aging.
the approach we followed uses stochastic transitions and it is showed in fig.
<dig>  where we represent a case that shows how th cells aging can be characterized.
here allow the binding of the variable th with any token with a product colorset th, we check that its life counter is greater than  <dig> with the guard expression  and, if such holds, we consume the input token and produce, by using the arc expression 1’, a new token that has the same state , the same position , the same receptor  of the input token, but with a life counter decreased by one .fig.
<dig> example of the use of a stochastic transition to model entities aging.
the guard on the stochastic transition randomly selects only tokens that represent th cells whose life is greater than  <dig> , and gives back to the starting place a new token that is identical to the previously selected one, except for the fact that life has been decremented by  <dig> )
interaction among two or three different entities represents another example of entity internal state changes.
as a result of the interaction of two or more other entities, cells may change their internal state and behavior, molecules like antigens may disappear from the system, and new entities can be introduced.
an example is showed in fig.
<dig>  where we give a possible representation of the interaction between primed th cells and presenting b cells.
both cells must be in the right internal state to interact; moreover receptor matching is needed to make the interactions possible.
as a result of a successful interaction  we will obtain duplication of th cells, and asymmetrical division/differentiation of b cells that will produce a memory b cell and/or a new plasma b cell.fig.
<dig> example of the use of a stochastic transition to model interaction among entities.
the guard on the stochastic transition randomly selects only tokens that represent th cells that are in the active state  and b cells that are in the presenting state .
moreover the activation of the transition is achieved only if both the cells have the same receptor specificity .
as a result of the activation of the transition we will obtain duplication of th cells by giving back  <dig> identical th cells to the th place , a memory b cell with the same position and specificity, but with renewed life )) and, as a result of a asymetrical division, a new plasma cell  with same receptor and characteristics of the progenitor b cells )
in particular, any token coming from the th and b places will be a potential candidate to be selected for the binding of the th and b variables, respectively.
however, the transition will be enabled only if the guard condition  that imposes that: 1) the token representing the b cell is in the “presenting” state; 2) the token representing the t helper cell is in the “active” state and, 3) there is matching within the involved receptors  is satisfied.
as result of the firing of the transition: 1) the input tokens will be consumed; 2) two newborn t helper cells with the same state, position and receptor of the input t helper cell, but with an increased life will be produced ) into the th place; 3) a newborn plasma b cell token with the same receptor and position of its progenitor, but with life increased to the maximum value will be introduced into the pb  place ); 4) a newborn memory b cell token with the same position and receptor of the input b cell, but with life increased to the maximum value will be produced into the b place ).
one of the most important features of the adaptive immune system response is given by the ability to recognize specific targets and characteristics of potential pathogens.
this goal is achieved by is cells using cell receptors that are able to recognize and/or capture specific amino acidic sequences that are present on the surface of infected cells and bacteria.
such a recognition process can be introduced directly inside transition guards and its definition depends on the receptor color definition itself.
in our example, since we only used a limited set of integer numbers to represent receptors, we supposed that a successful recognition could be achieved only if entities have the same integer number inside the receptor color.
of course, more complex models of both receptors and receptor recognition functions, like those based on binary strings and hamming distance implemented in the abm models presented in  <cit> , can be achieved by using specific function definitions based on look-up tables and more complex guards.
entity position and movement represents another aspect that can be taken into account for modeling is.
in some cases the presence of different compartments or an accurate spatial description can be useful to represent, for example, migration of is cells towards different compartments , tumor growth, and the spread of viral infection.
spatial movement can be achieved as described earlier for entity interactions by using stochastic transitions and guards.
in fig.
<dig> we show an example on how to implement random walk on two dimensions over a square lattice by using two colorsets x and y for position.
of course more complex spatial descriptions, like the one presented in  <cit>  for describing reaction diffusion systems, can be used.fig.
<dig> example of the use of a stochastic transitions to represent 2d random walk over a square grid with toroidal geometry.
in this example we suppose that the place that represents th entities contains tokens defined according the following product colorset: , with x and y defined as integer colorsets.
random 2d movement is achieved by using stochastic transitions.
the x_inc and x_dec transitions are used to model movement under the x direction, whereas the y_inc and y_dec represent the movement on the y direction.
the predecessor  and successor  colorset functions are here used to decrement or increment the position indices
we note here that in our proof-of-concept model, even if we predisposed the presence of a position colorset, we did not use it since we did not include any spatial description.
it must be highlighted that if space is included inside the model, all transitions guards that are used to model interactions must be revised accordingly to include such an aspect .
delayed events may represent another event that may be critical to represent.
for some scenarios it would be of major importance to represent events that require some time to complete, such as antigen internalization and processing by antigen presenting cells like macrophages.
such events can be represented making use of delayed transitions that allow the definition of one or more delay times for the activation of the considered transition.
biology implemented into the model
in order to demonstrate the feasibility of the modeling approach, we developed a model of the adaptive immune response whose goal is to reproduce two important is features: specificity and memory.
these two well-known features allow the immune system to recognize and attack specifically a given pathogen thanks to the use o specific cellular receptors , and to have an enhanced secondary response against a previously encountered pathogen .
in this paper we focused on the humoral immune response, that is driven by b lymphocytes and whose final outcome is represented by the production of antibodies.
we used snoopy  <cit>  to develop the stochastic colored pn according to the approach described earlier.
the model is shown in fig.
<dig> and is provided online .fig.
<dig> model of the humoral immune system response.
white circles represent places ; white boxes represent stochastic transitions; grayed boxes that represent scheduled transitions are used to introduce bacteria population at given times; the two black boxes that represent delayed transitions are used here to mimic the time needed by macrophages  and b cells to internalize, process and expose antigenic peptides
the biology implemented into the model can be reassumed as follows.
let us have a given pathogen inside the host i.e., an extracellular bacteria, that is able to replicate.
of course bacteria, like the other cellular entities, are exposed to aging and death.
in order to test the presence of the memory and specificity features of the is response, we supposed two subsequent infection spreads of two different bacterial papulations  into the system at given times using scheduled transitions.
bacteria can be recognized a-specifically by antigen presenting cells like macrophages  through the use of toll-like receptors or also specifically by resting b cells.
as a result of such interactions mp that will phagocyte bacteria will be activated and will proceed to process and then expose bacterial peptides complexed with major histocompatibility complex class ii peptides.
on the other hand, resting b cells that are able to specifically recognize bacterial determinants will be activated and release low affinity antibodies .
such antibodies will bind with low affinity to bacterial surfaces.
however, such a response may be not sufficient to eradicate bacteria.
to this end other mechanisms of the is response are involved.
activated mp will present exposed peptides to specific resting t helper cells  that will be primed; primed b cells that will encounter again bacterial antigens will be further stimulated and will become ready to present bacterial antigens to primed th cells.
when such interaction holds, th will be duplicate, and b cells will asymmetrically differentiate into plasma b cells  and/or memory b cells.
pb will produce specific antibodies  that will bind with higher affinity to bacterial surfaces, having as a result the death of bacteria .
memory b cells will indeed live longer and will be ready to directly differentiate into plasma cells if future encounters with the same bacteria happen.
all interactions can be examined in depth by analyzing the snoopy pn model we provided as supplementary material .
RESULTS
immune responses are specific for distinct antigens that are specifically recognized by individual lymphocytes.
the parts that the competent immune cells recognize are called determinants or epitopes.
this fine specificity exists because individual lymphocytes express membrane receptors that are able to distinguish elusive differences in structure between different epitopes.
the repertoire represents the total number of antigenic specificities of the lymphocytes in an individual.
it is extremely large and it has been estimated that the immune system of an individual can discriminate between  <dig> and  <dig> distinct antigenic determinants.
exposure of the immune system to a foreign antigen enhances its capacity to respond again to that antigen.
responses to second and subsequent exposures to the same antigen, called secondary immune responses, are usually more rapid, larger, and often qualitatively different from the first, or primary, immune response to that antigen .
immunologic memory occurs because each exposure to an antigen generates long-lived memory cells specific for the antigen, which are more numerous than the naive cells specific for the antigen that exist before antigen exposure.
in addition, these memory cells have special characteristics that make them more efficient at responding to and eliminating the antigen than are naive lymphocytes that have not previously been exposed to the antigen.
figure  <dig> sketches memory b lymphocytes that produce antibodies that bind antigens with higher affinities than do antibodies produced in primary immune responses.
in such an example antigen x is inoculated at the beginning.
immune system competent cells recognize the epitopes of the antigen x and mount the primary immune response, producing specific antibodies against antigen x. moreover, memory b cells are also produced.
later, another administration is performed.
this time, two different antigens are administered, antigen x and y. the presence of immune memory is highlighted by a larger response.
the antigen x is rapidly cleared while the response to the antigen y is slower, as this is the first time that immune competent cells face with its epitopes.fig.
<dig> typical immune system primary and secondary response.
different antigens  induce the stimulation of b cells with different specificity, and thus entitle the production of different antibodies .
the secondary response to the antigen x is more rapid and stronger than the primary response  and is different from the primary response to the antigen y .
antibody levels decline with time after each immunization
these features of adaptive immunity are necessary if the immune system is to perform its normal function of host defense.
specificity and memory enable the immune system to mount heightened responses to persistent or recurring exposure to the same antigen and thus to combat infections that are prolonged or occur repeatedly.
our goal is to demonstrate how the proposed modeling approach based on cpns can be useful for modeling the is response.
to this end, we developed a proof-of-concept model that is able to qualitatively reproduce the typical primary and secondary is response that is well known in the literature  <cit> , highlighting the presence of memory and specificity.
while the typical behavior of the immune system response described above is well known and represents one of the hallmarks of immunology, to our best knowledge there exists no specifically validated gold-standard clinical dataset that can be used to quantitatively compare, in a numerical sense and in a precise way, the outcomes of an in vivo/in silico experiment about the general immune system response to a given pathogen.
this is due to the fact that the immune system responses are influenced by multiple factors  that can entitle very high fluctuations even among different individuals of the same experimental group.
as a result of this, mathematical and computational models are usually tailored to the specific experiment they want to reproduce, and experimental data is usually used to fine-tune the model behavior.
in this case, since our goal is to demonstrate how high-level pns can be successfully applied to qualitatively reproduce specific immune system features rather than numerically reproduce a given experimental setup, we will show how memory and specificity features are coherently reproduced by qualitatively comparing pn based simulation results with those obtained using the uiss simulator, an agent based simulation framework that has been developed for years and specialized to model and predict the is response to different pathologies such as mammary carcinoma  <cit> , lung mestastases  <cit> , melanoma  <cit> , atherosclerosis  <cit> , and to predict and to optimize candidate treatments and novel adjuvants  <cit> .
the time step of all the simulations is Δt = 8 h, as in the abm models presented in  <cit> .
the total time of all simulations is  <dig> .
we recall here that the typical adaptive immune system response requires from  <dig> to 4 weeks.
we used the scheduled transitions to simulate the appearing of bacteria at different times.
since the initial marking of the pn has no tokens, before introducing any bacteria token, we waited a minimal number of time-steps  to give time to let the pn reach an equilibrium “homeostatic” level for the mp, th and b places .
all the experiments we have done have been repeated  <dig> times and mean values have been taken into account.
as a first experiment we challenged the model to test the induction of memory and compared these first results with those obtained with the uiss framework.
to this end we simulated the appearing of a given bacteria ) at two different time-steps.
we recall here that we did not include the state colorset in the bacteria definition.
the first appearing of the bacteria  is done to induce a primary response, whereas the second one  is done to test the presence of a stronger secondary response.
as similar setup  has been used for the uiss framework also.
results are shown in figs.
<dig> and  <dig>  looking at these figures, it is possible to see how the two modeling approaches show, from a qualitative point of view, show a similar behavior in reproducing the immune system response.
in particular, it appears clear that at the time of the second bacteria administration, the total quantity of igg is in both the models  and ) much higher than the quantity measured at the first injection.
another important aspect is represented by the fact that, as a consequence of the stronger immune response and of a higher peak of igg, the bacteria are killed more rapidly.
of course some aspects like the total number of antibodies needed and the  time required to clear the second bacterial infection are probably better reproduced by the uiss framework that, after years of work, contains a lot more of biological detail  and a spatial description that, at the present time, has not been enabled in the pn model.fig.
<dig> detailed dynamics of specific immune system entities in the pn model  and the uiss framework .
at the time of the second bacteria injection , the total level of igg  is much higher in respect to the level measured at the time of the first injection.
as a result of this, bacteria  are cleared more rapidly at the time of the second injection.
both igm  and igg  levels obtained with the pn model  are in good agreement  with the results obtained by using the uiss framework
fig.
<dig> detailed dynamics of specific memory b cells in the pn model and in the uiss computational framework.
memory b cells levels obtained with the pn model are in good agreement  with the results obtained by using the uiss framework.
the drop in the number of memory b cells in the pn model just after the second injection is caused by the fact that, after the second injection of bacteria, all memory cells are temporarily enrolled into the interaction with bacterial antigens, and rapidly move inside the relative transition .
this leaves for a short time the relative place without memory cells.
as a result of this, since we actually plot only tokens that are in places, for a very short time the relative curve goes down
once that we have checked that the pn model showed an is response that is qualitatively in line with the uiss framework, we tested the presence of the specificity feature by injecting two different bacterial populations .
the first injection of  <dig> tokens with the binding  is done at time  <dig>  and then a second injection of both the bacteria populations  and  <dig> tokens with binding , respectively) is done at time  <dig>  in other words, we reproduced the typical scenario of the immune system response that we have already presented in fig.
<dig>  in fig.
<dig> we show the results obtained with the pn model.
as expected, both specific plasma b cells and igg are higher for the first bacteria population ) at the time of the second injection.
furthermore, the response to the second bacteria population ) is comparable to the first population at the first time injection.
this highlights how the first injection elicits the formation of memory b lymphocytes that are specific for the first bacteria population only, and thus elicit a stronger immune response for that specific bacteria only.fig.
<dig> detailed dynamics of plasma cells and igg.
both specific plasma b cells ) and igg ) are higher for the first bacteria population  at the time of the second injection.
furthermore, the response to the second bacteria population  is comparable to the response obtained for the first population  at the first time injection.
this highlights how the first injection elicits the formation of specific memory b lymphocytes for the first bacteria population only
CONCLUSIONS
traditional “low-level” pn represent a strong modelling technique for the simulation of networks, concurrent systems and distributed systems in general.
even if pn have been initially used by computer science and engineering theoreticians only, thanks to the introduction of many extensions of the pn framework , pn have also been successfully applied in many other fields, e.g., for the modeling of biochemical pathways.
their initial application to the modeling of the immune system response has been taken into account.
however, such an approach did not make inroads, probably because it has not brought any real additional advantage over traditional top-down approaches.
for example, the lack of granularity in the description of the biological phenomena may have represented an important missing feature.
colored pn are a pn extension which may allow to overcome this limit.
to demonstrate this, we proposed a methodology to realize detailed models of the immune system response based on  cpn, and, to validate our methodology, we presented a model of the humoral adaptive is response against an extracellular bacteria that is able to reproduce well-known important is features such as memory and specificity.
the proposed methodology has some advantages over classical methodologies.
it allows to come close to the granularity allowed by agent based models, by letting the modeler to describe in depth, for example, internal cell states, receptors and complex interactions.
at the same time cpn benefit of many types of formal analysis that may permit to analyze some fundamental properties of the biological network under study such as reachability, liveness, boundedness, place invariants, circuits, t-invariants and so on  <cit> , making them more similar to differential equations models in this sense.
furthermore, cpn represent a graphical modeling technique that does not require a strong knowledge of mathematics, like in differential equation based models, or good programming skills, like for abm models.
thus, cpn based models are easy to design, pre-existing software like snoopy or cpn tools can be used, and life-scientists that already demonstrated able to use pn tools in general for the developing and analysis of many biochemical signaling pathway models may use cpn with little or no effort at all.
indeed, an additional benefit of the use of cpn is represented by the possibility to realize multi scale models, but with the advantage of using the same modeling technique and software to combine, using hierarchical petri nets, intra cellular models  and system models .

 the immune system is a highly distributed system composed by cells and molecules whose scope is to protect living organisms against diseases. the immune system has various defense mechanisms whose complexity is usually correlated to the host organism. in mammals, among many defense mechanisms such as innate immune response arm, the adaptive is response, driven by lymphocytes, represents a fundamental component. the goal is of the adaptive is response is to give the best reaction against most dangerous viruses and bacteria. specificity, memory and discrimination between self and non-self represent the most important peculiarities of the adaptive is response. in particular, the specificity feature is the capability, shared only by lymphocytes and antibodies, to recognize one specific epitope of foreign pathogens. specific immune system lymphocytes can then recognize billions of different antigens, but for each determinant, a specific lymphocyte will be stimulated. there are as many stimulated lymphocytes as determinants forming the antigen. indeed, the memory feature represents the capability of the immune system to remember previously encountered antigens. when an antigen appears for the first time, the immune system produces a first response against it. this stimulates the production of long-lived memory lymphocytes that are specific for each of the antigen epitopes. when the same antigen comes again into contact with the immune system, the memory lymphocytes are stimulated in order to produce a faster and more effective immune system response against it. building a model of the immune system may lead to multiple benefits, since researchers can use immune system simulation in vaccine research, or to test biological hypotheses about diseases and infection dynamics. many mathematical and computational techniques have been applied to model the immune system response; most of them bear on two main approaches: the top-down approach and the bottom-up approach. if the former mainly focuses on the behavior of a system as a whole, the latter analyzes the system dynamics at a lower level and the emergence of global properties is obtained as the sum of local behaviors. models based on the top-down approach usually make use of differential equations like ordinary differential equations. these models are all population-based, and the spatiality and topology which both depend on individual interactions are, in general, ignored. for very simple models it is possible to use some powerful mathematical techniques to extract some important system properties like analytical solutions, steady states and asymptotic behavior. however, as the biological description introduced into the model grows, they become very complex and, in some cases, almost intractable. indeed, the application of partial differential equations, stochastic differential equations and delayed differential equations can be introduced to obtain some spatial description, stochastic effects and temporal delays. however, this usually entitles instability problems, a higher computational effort and numerical simulations only. cellular automata and agent based models represent the most famous bottom-up approaches. since entities are followed individually, they allow a much more rich description of the biological background, and individual cells properties like cell receptors, cell internal states and individual life-times and behaviors are not a problem. furthermore, spatial description, stochastic effects and delays are intrinsic characteristics of these kind of models. however, such techniques lack of tools and mathematical techniques for analytical studies and they usually require huge computational resources for large scale simulations. petri nets represent a graphical and mathematical modeling technique for the description of distributed systems presented in by carl adam petri. during these last decades, petri net modeling formalism has been extended in order to manage, among others, stochastic effects, continuous quantities and scheduled events. among pn extensions there is the colored petri nets formalism. petri nets and extensions have been successfully used in bioinformatics for the modeling of intra-cellular biochemical and metabolic pathways. park et al. showed how it is possible to apply a particular extension of pn for the modeling of the immune system response. the approach had several limitations; among them, it did not allow to discern between cells featuring different internal states, receptors, or a different position inside the simulation space. consequently, such an approach did not bring any real additional advantage over up-down traditional techniques such as ode systems, and neither in respect to abms. here we propose a novel methodological approach that uses “high level” petri nets, and in particular colored petri nets to model cell specific properties such as cell internal states and receptors. to this end we present a proof-of-concept model of the humoral part of the is capable of reproducing the most important features of the adaptive response, and in particular the memory and specificity features. the approach holds somewhere in the middle of top-down and bottom-up approaches, since it allows to describe cell specific features and behaviors, but with the possibility of using many standard pn analysis tools that permit to study system properties and features from a qualitative point of view. all of this using a convenient graphical language that is already known to life scientists that, as already stated, are actually able to use pn for the modeling of signaling pathways. the paper starts with an introduction to petri nets and the colored extension used in this paper. details of the modeling approach are provided in methods, followed by results and discussion where we present the application of the approach to a real scenario, ending with our conclusions and final remarks. methods introduction to classical “low-level” petri nets pns is a set of particular graphical symbols called: places, transitions, arcs and tokens. the places are generally drawn as ellipses or circles and represent states and/or conditions of the system modeled; each place may be featured by a name, which has not a formal meaning but increases understanding of a petri net. boxes or rectangular bars represent transitions which model events that may occur in the system. places usually contain a discrete number of elements called tokens, drawn by black dots. arcs may connect places to transitions and vice-versa, but they cannot connect transitions to transitions or places to places. an arc connects a single place to a single transition; furthermore, an arc may connect a transition to a single place. a transition may have several input and output places. the tokens in the input places are fundamental for the firing of the transition; the enabling of a transition models the occurrence of a specific event. firing of a transition may require the presence of a different number of tokens for each input place. generally, the required number of tokens for each input place is given by the multiplicity of each arc connecting input places with the transition. when a transition fires, it consumes the required number of tokens from the input places, and creates tokens in its output places. the number of tokens put in each output place depends again by the multiplicity associated to each arc. a firing rule is atomic, in the sense that it represents a single non-interruptible action. unless an execution policy is defined, when multiple transitions can be enabled at the same time anyone of them may fire. introduction to “high-level” colored petri nets a colored petri net extend the definition of petri nets as it is a set of places, transitions, arcs and color sets. each place has an associated type determining the kind of data that the place may contain. like pns, during the execution of a cpn each place will contain a varying number of tokens. each of these tokens carries a data value that belongs to the type associated with the place. for historical reasons, token values are referred as token colors and data types are referred as color sets. this is a metaphoric picture where the tokens of “high-level” cpn are considered to be distinguishable from each other and hence “colored” – in contrast to ordinary “low-level” petri nets which have “black” indistinguishable tokens. the types of a cpn can be arbitrarily complex e.g., a record where one field is a real, another a text string and a third a list of integers. hence, it is much more adequate to imagine a continuum of colors instead of a few discrete color values. in cpns, there are different types of expressions. an expression is mainly made up by variables. it is not only associated with a particular color set but also written in terms of a predefined syntax. in the following, we denote by exp a set of expressions that comply with a predefined syntax. an expression may be associated to a transition; in this case, it is called guard expression. an arc may have an expression associated, too; it is called arc expressions. on the basis of what said until now, the formal definition of colored petri nets will be given in the following. it uses the concept of multiset, which can be defined as a set in which there can be several occurrences for the same element. the number of occurrences of an element is called coefficient or multiplicity. for example, an infinite number of multisets exist which contain elements a and b, varying only by multiplicity, such as: {a, b}, {a, a, b}, {a, a, a, b, b, b}. in the following, given a generic set e, the collection of all the multisets over e is denoted by e ms. definition: a colored petri net is a tuple n =  where:s is a finite non-empty set of types, called color sets. p is a finite non-empty set of places. t is a finite non-empty set of transitions. a is a finite set of direct arcs c is a color function that assigns to each place p∈p a color set c∈s. c: p → s g is the set of g functions that assigns to each transition t∈t a guard expression of the boolean type. g: t → exp f is the set of f function that assigns to each arc a∈a an arc expression of a multiset type cms, where p is the place connected to the arc a. f: a → exp i is an initialization function that assigns to each place p∈p an expression of a multiset type cms. i: p → exp according to the tokens distribution in the set p∈p, the cpn may assume different states; each state is called marking and is given by a particular token distribution in each of the places p∈p. the initial marking is produced by the application of the initialization function i. in the following the current token distribution in a particular place will be indicated by m, which is a function of a multiset type cms. given a particular a transition t ∈ t, let us consider the set of variables present in the guard of t and in the arc expressions of arcs connected to t; this set will be denoted var. before the guard and arc expressions relevant to the transition t ∈ t may be evaluated, the variables must get assigned values; this process is called binding. a binding b of a transition t ∈ t is a function that maps each variable v ∈ var onto a value b that is of the same type as the variable. the set of all bindings for a transition t is denoted b. a binding element is a pair, with t ∈ t and b ∈ b. the behavior of the cpn is based on a firing rule consisting of a precondition and the effect of the occurrence of a single transition. whether or not a transition can fire depends on the marking of its input places and the arc expression on the input arcs. a transition is enabled and is allowed to fire, if all the input places are sufficiently marked; this occurs when the binding of the variables that appear in the arc expressions of the arcs from inputs places to transition, evaluates to a multiset of token colors ms) that is present on the corresponding input place. in this case, the guard of the transition should evaluate to true for the giving binding. if a transition has no input places, it is always enabled. in a more formal way, it is possible to say that the guard of the transition g ∈ g evaluates to true if m ≥ f, ∀p ∈ •t, where f ∈ f and •t is the set of input places for transition t. due to the definition of set f, f denotes the value ms) of the arc expression relevant to the arc going from p to t, for the giving binding b. if the guard of the transition evaluates to true, the firing of a transition occurs. in this case, a multiset of colored tokens are taken from each input place and added to output places in accordance with the values given by the arc expressions relevant to the arcs connecting the transition to output places. this means that the firing process lead to a change in the current marking of the cpn. in particular, the current marking m becomes the new marking m’, defined by m’ = m + f’-f, where f’ ∈ f is evaluated for ∀p ∈ t• and denotes the values of the arc expressions relevant to the arcs going from p to t, for the giving binding b. a binding element may occur concurrently to other binding elements – iff there are so many tokens that each binding element can get its “own share”. the colored generalized stochastic petri net framework extension different cpn modeling paradigms have been defined during years: colored qualitative petri nets, colored stochastic petri nets, and colored continuous petri nets. qpn c does not uses tokens but it is only used for qualitative analysis of systems; spn c takes into account discrete tokens and stochastic transitions, and cpn c uses continuous quantities rather than discrete tokens and continuous transition rates instead of discrete transitions. qpn c is an abstraction of spn c and cpn c, while spn c can reproduce by approximation cpn c and vice-versa. from our point of view, the most interesting approach is represented by spn c. in brief, this colored pn extension uses a discrete number of tokens on its places like in classical pns. moreover, stochastic events and transitions that are fundamental for simulating, for example, different individuals, are flawlessly reproduced through the use of stochastic transitions which fire after a probabilistic delay determined by a random variable defined according an exponential probability distribution. for our purposes it is useful to consider a further extension of the spn c framework, namely colored generalized stochastic petri nets, by considering, besides of stochastic transitions, also some kinds of deterministic transitions that we can identify as immediate, delayed, and scheduled transitions. all transitions become enabled if all the preplaces are sufficiently marked. when a stochastic transition is enabled, a given time must be wait before the firing occurs. this waiting time that determines the firing delay of the transition is given by a random variable x t that is distributed exponentially with the following probability density function: fxtτ=λtm⋅e−λtm⋅τ,τ≥ where λt can be defined as an arbitrary mathematical function which depends on the marking m of the preplaces at time t. it should be noted that, even if there is a stochastic time delay before the firing of the transition, the firing itself does not consume any time. deterministic transitions have a deterministic firing delay specified by an integer number. the delay count starts just after the transition becomes enabled. however, it must be said that during this waiting time it can happen that the transition loses its enabled state. immediate transitions can be seen as a special case of deterministic transitions with a delay time set to in case of conflict between an immediate transition and any other kind of transition, the former will get firing priority. also scheduled transitions can be seen as a special case of deterministic transitions. the firing is deterministic and it occurs at a previously defined absolute time of the simulation, obviously only if the transition is enabled at that time. interested readers can find further details about spn c. advantages of using colored petri nets petri nets represent a graphical modeling tool that allows to describe in a simple and clear, but yet formally correct and powerful manner, any kind of process. the biggest problem of classical low-level petri nets is given by the fact that they usually not scale. as the process we want to describe grows, low-level petri nets tend to grow quickly, and the designing, drawing, understanding and managing of the net becomes even more and more difficult, thus increasing the developing time and the risk of introducing modeling errors. the introduction of “colors” allows to discern among different tokens and to follow the individually, in a similar way as it is done by ca and abm approaches. this represents the most important feature of colored petri nets, and their greatest advantage in respect to low-level pns, but also in respect to all the differential equation based models. just to give to the reader an idea of the practical advantage of cpns, the proof-of-concept model that we developed according to the approach described in the next section required only places, transitions and arcs. if a standard pn approach would have been used, it would have required approximately different places, transitions and arcs. this means that the development of the model itself would have been almost unfeasible. a similar scenario holds, for example, if a standard ode approach is used. the biggest advantage of cpns in respect to classical bottom up approaches is given by the availability of a large number of formal analysis methods by which the important properties of the net can be extrapolated, proved and analyzed. liveness, boundedness, place invariants, circuits, t-invariants and so on, allow to study the net properties from a qualitative point of view, making the pn approach more similar to a differential equation based approach in this sense. furthermore, over classical simulation, it is possible to have interactive simulations with results drawn directly and in real time on the cpn diagram, and it is possible to act during these simulations by, for example, introducing and modifying new and existing tokens. this kind of interactive simulation allows to debug the model behavior while it is developed, just like in a complete programming language environment. indeed, novel simulation techniques that should speedup and potentially allow large scale simulations as, for example, simulations that use the “symbolic marking” of the net or “fluidification” techniques have now been presented and their implementation in present pn software is in due course. finally, cpns have a graphical representation. this graphical representation is very simple, intuitive, and very appealing, even for people who are not very familiar with the formal details of cpn. with only three kinds of graphical items it is possible to mostly describe any kind of process and to obtain, for free, a very clear and qualitative conceptual model of process itself. arc and transitions guards can include very simple first-order logic rules, as well as more complex functions that can be defined using look-up tables. thus, nor a strong knowledge of mathematics, neither a good knowledge of any programming language is required, and this may represent a consistent leap ahead in respect to classical approaches, since it could allow to life-scientists to directly use cpns with little or no effort at all. how to use cpn to model the immune system function: a novel methodological approach to develop the adaptive is response model, we used colored generalized stochastic petri nets, in which, as already described, the most common transition rules are defined as stochastic processes. we used the snoopy software for developing a proof-of-concept model of the humoral response to an extracellular bacterium. the approach we present can be reassumed as follows:use places to identify cellular or molecular types, or subpopulations types use transitions to model events, such as cellular and molecular state changes, entity interactions, duplication, differentiation, aging and death. use tokens to represent entity instances, like cells and molecules. use colors to define specific properties of cellular and molecular entities. in particular, we defined two distinct colorset types that include the following colours: state, position, lifetime and receptor. the state colorset is used to define internal cell states that determine cells behavior and it is actively used only for cellular types. in the present implementation its possible values are: “resting; active; presenting; memory; plasma”. the life colorset is used to define life-time of cellular and molecular entities entities. the receptor colorset can be used for many functions and its meaning depends on the cellular/molecular entity and transition rule we are dealing with. it can represent:the cellular receptor that is used for the recognition of specific targets by b or t cells; the presented peptide sequence that is presented by antigen presenting cells like mp; the epitope sequence presented bacteria or antigens in general; the igg antibody specific binding structure. the position colorset can be used to model position of cellular entities inside the space or to represent different compartments, i.e. viral infection place and lymph node. by using this approach and colored petri nets it will be possible to describe in a very general way, with just one transition, the behavior of all the entities that behave similarly, forgetting about the properties that are not relevant for the behavior. in other words it will be possible to say, for example, with just only one transition: “all the naive b cells that encounter in the same position and recognize with sufficient affinity a specific pathogen will be activated, no matter what is the valule of their real exact position, their value of their real exact receptor or peptide sequence, and how much remaining life they have”. with this approach we will gain a more grained description that introduces, for example, cell receptors and allows to follow entities individually, but at the same time we will obtain a very easy but yet powerful and flexible way to define the general rules that drive the dynamics of the model. of course, the definition of the used colors can be adapted to the problem we are dealing with i.e., the position colorset can be removed if we are not interested in modeling space and/or different compartments. on the contrary, we can consider two or tree different position colors to describe two or three dimensional environments. in our proof-of-concept model we have taken into account different places that are used to represent bacteria, macrophages, t helper cells, b cells, plasma cells, immunoglobulins class m and immunoglobulins class g. as already stated, we defined two distinct product colorsets types for all tokens, one type that is mainly used for cellular entities and one type that is used for molecular entities. this first product colorset is composed by four colors: state, position, lifetime and receptor, and it is used by th, b, pb, mp; the second one is composed by three colors: position, lifetime, and receptor and it is used for igm, igg, and bacteria. we note here that, even if bacteria are cells and not molecules, we used the second colorset definition for such entities since we were not interested in describing bacteria internal states. it is worth nothing that in the present model implementation all cells and molecules have distinct product colorsets in order to improve readability. however, it must be noted all the cellular entities have the same product colorset structure, and this also holds among molecular types. therefore, it would be possible to define an unique product colorset for all cells and another one for all molecules, since the distinction among different cell types or different molecular entities depends on the context. in fact, since places represent entity types, two identical tokens may represent two different kinds of cells if they hold in different places. the transition rules that are based on well-known biological rules define the correct state changes of the system and do not allow the erroneous migration of tokens towards places that represent different entity types. it will be assumed that for any arc expression we can either have a variable with of a given product colorset type that will bind to any token with the same colorset, or a more complex expression, always based on the colorset definition that will select only tokens that will satisfy some given particular requirements. for example the arc expression 2’, with p representing a variable whose type is position and l a variable whose type is life, when placed on an input arc will only accept two tokens with subset colors state=“active” and receptor=“4”. p and l will bind to any position and life colors of the incoming tokens. in other words, any token representing a cell in the active state and receptor will be accepted by the arc expression, no matter what life and position it will have. the same arc expression when defined on an output arc will produce two tokens with colors state=“active”, receptor=“4”, and position and life colors defined according to the p and l values. if p and l binding already occurred, for example on the input arcs of the transition, their value will be already defined and assigned on the tokens that will be produced; if no binding of p and l occurred on the input arc expressions, their value will be chosen at random inside the set of the allowed position and life colorset values. this example takes into account tokens which represent cellular entities, however, if we remove the state color set from the product color set definition, it is also valid for molecular entities. as already stated, the evolution of the system is driven by transition rules that are defined upon well-known biological rules. generally, transitions can be used to delineate all the biological scenarios that are needed to describe the evolution of the system. among such scenarios we recall the appearing of new entities and death of old entities, cellular and molecular behavior changes that include changes in the internal cell state, aging, processing and presentation of antigenic sequences and movement. such changes can be the result of transitions that involve one or more entities. of course it is not possible to describe in detail all the possible biological scenarios that may arise, so we sketch out more explicitly the methodology by showing the most relevant examples. introduction and death of new entities the introduction and of new entities can be modeled as presented in fig.  usually, new entities can be introduced with stochastic transitions with an initial pre-defined color. for example, as shown in fig.  th cells are introduced as “resting” cells, with an initial life-time set to both the receptor and the position colors are not assigned, and this translates into the generation of new cells with random position and receptor, but with a predefined life counter and initial internal state. as it also possible to note from fig.  such kind of transitions does not need any incoming token from a given place to be enabled.fig. example of introduction of new entities into the system. new th cells are introduced at random times into the system as resting cells, with a random position, life =  and a random receptor for some scenarios it is requested to model the introduction of new entities at a given time, for example to reproduce the injection of a given treatment or compound. if such situation holds, it is possible to substitute the simple transition with a scheduled transition that allows the definition of an initial time, period, and final time of activation of the transition. death of entities will occur for all cells and/or molecules that have ended their life cycle. this can be achieved by using transitions with transition guards on the value of the component element that refers to the life counter, as presented in the example shown in fig.  here, any token can be potentially selected, in fact any token coming from the th place can bind to the variable th, but the transition will fire only if the guard condition on the transition is satisfied, that is if the third component of th, representing the life subset color, is if such holds, the input token will be consumed by the transition, and no output tokens will be produced.fig. example of the use of a stochastic transition to model the dead of entities that have completed their life cycle. th cells that have completed their life cycle will be removed from the system. the transition is enabled if and only if the guard condition on the life subset color holds, thus only if the binding on the th variable is done on a token whose third component element of the product colorset th is equal to zero internal state changes with the preposition “internal state changes”, we define all the events that modify the binding of the components of the product colorsets we have previously defined. such changes may hold from an immunological point of view for multiple reasons. a first example of this scenario can be represented by entities aging. in this case the, life component of entities needs to be decreased to reflect aging. the approach we followed uses stochastic transitions and it is showed in fig.  where we represent a case that shows how th cells aging can be characterized. here allow the binding of the variable th with any token with a product colorset th, we check that its life counter is greater than with the guard expression and, if such holds, we consume the input token and produce, by using the arc expression 1’, a new token that has the same state, the same position, the same receptor of the input token, but with a life counter decreased by one.fig. example of the use of a stochastic transition to model entities aging. the guard on the stochastic transition randomly selects only tokens that represent th cells whose life is greater than, and gives back to the starting place a new token that is identical to the previously selected one, except for the fact that life has been decremented by ) interaction among two or three different entities represents another example of entity internal state changes. as a result of the interaction of two or more other entities, cells may change their internal state and behavior, molecules like antigens may disappear from the system, and new entities can be introduced. an example is showed in fig.  where we give a possible representation of the interaction between primed th cells and presenting b cells. both cells must be in the right internal state to interact; moreover receptor matching is needed to make the interactions possible. as a result of a successful interaction we will obtain duplication of th cells, and asymmetrical division/differentiation of b cells that will produce a memory b cell and/or a new plasma b cell.fig. example of the use of a stochastic transition to model interaction among entities. the guard on the stochastic transition randomly selects only tokens that represent th cells that are in the active state and b cells that are in the presenting state. moreover the activation of the transition is achieved only if both the cells have the same receptor specificity. as a result of the activation of the transition we will obtain duplication of th cells by giving back identical th cells to the th place, a memory b cell with the same position and specificity, but with renewed life )) and, as a result of a asymetrical division, a new plasma cell with same receptor and characteristics of the progenitor b cells ) in particular, any token coming from the th and b places will be a potential candidate to be selected for the binding of the th and b variables, respectively. however, the transition will be enabled only if the guard condition that imposes that: 1) the token representing the b cell is in the “presenting” state; 2) the token representing the t helper cell is in the “active” state and, 3) there is matching within the involved receptors is satisfied. as result of the firing of the transition: 1) the input tokens will be consumed; 2) two newborn t helper cells with the same state, position and receptor of the input t helper cell, but with an increased life will be produced ) into the th place; 3) a newborn plasma b cell token with the same receptor and position of its progenitor, but with life increased to the maximum value will be introduced into the pb place ); 4) a newborn memory b cell token with the same position and receptor of the input b cell, but with life increased to the maximum value will be produced into the b place ). one of the most important features of the adaptive immune system response is given by the ability to recognize specific targets and characteristics of potential pathogens. this goal is achieved by is cells using cell receptors that are able to recognize and/or capture specific amino acidic sequences that are present on the surface of infected cells and bacteria. such a recognition process can be introduced directly inside transition guards and its definition depends on the receptor color definition itself. in our example, since we only used a limited set of integer numbers to represent receptors, we supposed that a successful recognition could be achieved only if entities have the same integer number inside the receptor color. of course, more complex models of both receptors and receptor recognition functions, like those based on binary strings and hamming distance implemented in the abm models presented in, can be achieved by using specific function definitions based on look-up tables and more complex guards. entity position and movement represents another aspect that can be taken into account for modeling is. in some cases the presence of different compartments or an accurate spatial description can be useful to represent, for example, migration of is cells towards different compartments, tumor growth, and the spread of viral infection. spatial movement can be achieved as described earlier for entity interactions by using stochastic transitions and guards. in fig.  we show an example on how to implement random walk on two dimensions over a square lattice by using two colorsets x and y for position. of course more complex spatial descriptions, like the one presented in for describing reaction diffusion systems, can be used.fig. example of the use of a stochastic transitions to represent 2d random walk over a square grid with toroidal geometry. in this example we suppose that the place that represents th entities contains tokens defined according the following product colorset:, with x and y defined as integer colorsets. random 2d movement is achieved by using stochastic transitions. the x_inc and x_dec transitions are used to model movement under the x direction, whereas the y_inc and y_dec represent the movement on the y direction. the predecessor and successor colorset functions are here used to decrement or increment the position indices we note here that in our proof-of-concept model, even if we predisposed the presence of a position colorset, we did not use it since we did not include any spatial description. it must be highlighted that if space is included inside the model, all transitions guards that are used to model interactions must be revised accordingly to include such an aspect. delayed events may represent another event that may be critical to represent. for some scenarios it would be of major importance to represent events that require some time to complete, such as antigen internalization and processing by antigen presenting cells like macrophages. such events can be represented making use of delayed transitions that allow the definition of one or more delay times for the activation of the considered transition. biology implemented into the model in order to demonstrate the feasibility of the modeling approach, we developed a model of the adaptive immune response whose goal is to reproduce two important is features: specificity and memory. these two well-known features allow the immune system to recognize and attack specifically a given pathogen thanks to the use o specific cellular receptors, and to have an enhanced secondary response against a previously encountered pathogen. in this paper we focused on the humoral immune response, that is driven by b lymphocytes and whose final outcome is represented by the production of antibodies. we used snoopy to develop the stochastic colored pn according to the approach described earlier. the model is shown in fig.  and is provided online.fig. model of the humoral immune system response. white circles represent places; white boxes represent stochastic transitions; grayed boxes that represent scheduled transitions are used to introduce bacteria population at given times; the two black boxes that represent delayed transitions are used here to mimic the time needed by macrophages and b cells to internalize, process and expose antigenic peptides the biology implemented into the model can be reassumed as follows. let us have a given pathogen inside the host i.e., an extracellular bacteria, that is able to replicate. of course bacteria, like the other cellular entities, are exposed to aging and death. in order to test the presence of the memory and specificity features of the is response, we supposed two subsequent infection spreads of two different bacterial papulations into the system at given times using scheduled transitions. bacteria can be recognized a-specifically by antigen presenting cells like macrophages through the use of toll-like receptors or also specifically by resting b cells. as a result of such interactions mp that will phagocyte bacteria will be activated and will proceed to process and then expose bacterial peptides complexed with major histocompatibility complex class ii peptides. on the other hand, resting b cells that are able to specifically recognize bacterial determinants will be activated and release low affinity antibodies. such antibodies will bind with low affinity to bacterial surfaces. however, such a response may be not sufficient to eradicate bacteria. to this end other mechanisms of the is response are involved. activated mp will present exposed peptides to specific resting t helper cells that will be primed; primed b cells that will encounter again bacterial antigens will be further stimulated and will become ready to present bacterial antigens to primed th cells. when such interaction holds, th will be duplicate, and b cells will asymmetrically differentiate into plasma b cells and/or memory b cells. pb will produce specific antibodies that will bind with higher affinity to bacterial surfaces, having as a result the death of bacteria. memory b cells will indeed live longer and will be ready to directly differentiate into plasma cells if future encounters with the same bacteria happen. all interactions can be examined in depth by analyzing the snoopy pn model we provided as supplementary material.  immune responses are specific for distinct antigens that are specifically recognized by individual lymphocytes. the parts that the competent immune cells recognize are called determinants or epitopes. this fine specificity exists because individual lymphocytes express membrane receptors that are able to distinguish elusive differences in structure between different epitopes. the repertoire represents the total number of antigenic specificities of the lymphocytes in an individual. it is extremely large and it has been estimated that the immune system of an individual can discriminate between and distinct antigenic determinants. exposure of the immune system to a foreign antigen enhances its capacity to respond again to that antigen. responses to second and subsequent exposures to the same antigen, called secondary immune responses, are usually more rapid, larger, and often qualitatively different from the first, or primary, immune response to that antigen. immunologic memory occurs because each exposure to an antigen generates long-lived memory cells specific for the antigen, which are more numerous than the naive cells specific for the antigen that exist before antigen exposure. in addition, these memory cells have special characteristics that make them more efficient at responding to and eliminating the antigen than are naive lymphocytes that have not previously been exposed to the antigen. figure  sketches memory b lymphocytes that produce antibodies that bind antigens with higher affinities than do antibodies produced in primary immune responses. in such an example antigen x is inoculated at the beginning. immune system competent cells recognize the epitopes of the antigen x and mount the primary immune response, producing specific antibodies against antigen x. moreover, memory b cells are also produced. later, another administration is performed. this time, two different antigens are administered, antigen x and y. the presence of immune memory is highlighted by a larger response. the antigen x is rapidly cleared while the response to the antigen y is slower, as this is the first time that immune competent cells face with its epitopes.fig. typical immune system primary and secondary response. different antigens induce the stimulation of b cells with different specificity, and thus entitle the production of different antibodies. the secondary response to the antigen x is more rapid and stronger than the primary response and is different from the primary response to the antigen y. antibody levels decline with time after each immunization these features of adaptive immunity are necessary if the immune system is to perform its normal function of host defense. specificity and memory enable the immune system to mount heightened responses to persistent or recurring exposure to the same antigen and thus to combat infections that are prolonged or occur repeatedly. our goal is to demonstrate how the proposed modeling approach based on cpns can be useful for modeling the is response. to this end, we developed a proof-of-concept model that is able to qualitatively reproduce the typical primary and secondary is response that is well known in the literature, highlighting the presence of memory and specificity. while the typical behavior of the immune system response described above is well known and represents one of the hallmarks of immunology, to our best knowledge there exists no specifically validated gold-standard clinical dataset that can be used to quantitatively compare, in a numerical sense and in a precise way, the outcomes of an in vivo/in silico experiment about the general immune system response to a given pathogen. this is due to the fact that the immune system responses are influenced by multiple factors that can entitle very high fluctuations even among different individuals of the same experimental group. as a result of this, mathematical and computational models are usually tailored to the specific experiment they want to reproduce, and experimental data is usually used to fine-tune the model behavior. in this case, since our goal is to demonstrate how high-level pns can be successfully applied to qualitatively reproduce specific immune system features rather than numerically reproduce a given experimental setup, we will show how memory and specificity features are coherently reproduced by qualitatively comparing pn based simulation results with those obtained using the uiss simulator, an agent based simulation framework that has been developed for years and specialized to model and predict the is response to different pathologies such as mammary carcinoma, lung mestastases, melanoma, atherosclerosis, and to predict and to optimize candidate treatments and novel adjuvants. the time step of all the simulations is Δt = 8 h, as in the abm models presented in. the total time of all simulations is. we recall here that the typical adaptive immune system response requires from to 4 weeks. we used the scheduled transitions to simulate the appearing of bacteria at different times. since the initial marking of the pn has no tokens, before introducing any bacteria token, we waited a minimal number of time-steps to give time to let the pn reach an equilibrium “homeostatic” level for the mp, th and b places. all the experiments we have done have been repeated times and mean values have been taken into account. as a first experiment we challenged the model to test the induction of memory and compared these first results with those obtained with the uiss framework. to this end we simulated the appearing of a given bacteria ) at two different time-steps. we recall here that we did not include the state colorset in the bacteria definition. the first appearing of the bacteria is done to induce a primary response, whereas the second one is done to test the presence of a stronger secondary response. as similar setup has been used for the uiss framework also. results are shown in figs.  and looking at these figures, it is possible to see how the two modeling approaches show, from a qualitative point of view, show a similar behavior in reproducing the immune system response. in particular, it appears clear that at the time of the second bacteria administration, the total quantity of igg is in both the models and ) much higher than the quantity measured at the first injection. another important aspect is represented by the fact that, as a consequence of the stronger immune response and of a higher peak of igg, the bacteria are killed more rapidly. of course some aspects like the total number of antibodies needed and the time required to clear the second bacterial infection are probably better reproduced by the uiss framework that, after years of work, contains a lot more of biological detail and a spatial description that, at the present time, has not been enabled in the pn model.fig. detailed dynamics of specific immune system entities in the pn model and the uiss framework. at the time of the second bacteria injection, the total level of igg is much higher in respect to the level measured at the time of the first injection. as a result of this, bacteria are cleared more rapidly at the time of the second injection. both igm and igg levels obtained with the pn model are in good agreement with the results obtained by using the uiss framework fig. detailed dynamics of specific memory b cells in the pn model and in the uiss computational framework. memory b cells levels obtained with the pn model are in good agreement with the results obtained by using the uiss framework. the drop in the number of memory b cells in the pn model just after the second injection is caused by the fact that, after the second injection of bacteria, all memory cells are temporarily enrolled into the interaction with bacterial antigens, and rapidly move inside the relative transition. this leaves for a short time the relative place without memory cells. as a result of this, since we actually plot only tokens that are in places, for a very short time the relative curve goes down once that we have checked that the pn model showed an is response that is qualitatively in line with the uiss framework, we tested the presence of the specificity feature by injecting two different bacterial populations. the first injection of tokens with the binding is done at time and then a second injection of both the bacteria populations and tokens with binding, respectively) is done at time in other words, we reproduced the typical scenario of the immune system response that we have already presented in fig.  in fig.  we show the results obtained with the pn model. as expected, both specific plasma b cells and igg are higher for the first bacteria population ) at the time of the second injection. furthermore, the response to the second bacteria population ) is comparable to the first population at the first time injection. this highlights how the first injection elicits the formation of memory b lymphocytes that are specific for the first bacteria population only, and thus elicit a stronger immune response for that specific bacteria only.fig. detailed dynamics of plasma cells and igg. both specific plasma b cells ) and igg ) are higher for the first bacteria population at the time of the second injection. furthermore, the response to the second bacteria population is comparable to the response obtained for the first population at the first time injection. this highlights how the first injection elicits the formation of specific memory b lymphocytes for the first bacteria population only  traditional “low-level” pn represent a strong modelling technique for the simulation of networks, concurrent systems and distributed systems in general. even if pn have been initially used by computer science and engineering theoreticians only, thanks to the introduction of many extensions of the pn framework, pn have also been successfully applied in many other fields, e.g., for the modeling of biochemical pathways. their initial application to the modeling of the immune system response has been taken into account. however, such an approach did not make inroads, probably because it has not brought any real additional advantage over traditional top-down approaches. for example, the lack of granularity in the description of the biological phenomena may have represented an important missing feature. colored pn are a pn extension which may allow to overcome this limit. to demonstrate this, we proposed a methodology to realize detailed models of the immune system response based on cpn, and, to validate our methodology, we presented a model of the humoral adaptive is response against an extracellular bacteria that is able to reproduce well-known important is features such as memory and specificity. the proposed methodology has some advantages over classical methodologies. it allows to come close to the granularity allowed by agent based models, by letting the modeler to describe in depth, for example, internal cell states, receptors and complex interactions. at the same time cpn benefit of many types of formal analysis that may permit to analyze some fundamental properties of the biological network under study such as reachability, liveness, boundedness, place invariants, circuits, t-invariants and so on, making them more similar to differential equations models in this sense. furthermore, cpn represent a graphical modeling technique that does not require a strong knowledge of mathematics, like in differential equation based models, or good programming skills, like for abm models. thus, cpn based models are easy to design, pre-existing software like snoopy or cpn tools can be used, and life-scientists that already demonstrated able to use pn tools in general for the developing and analysis of many biochemical signaling pathway models may use cpn with little or no effort at all. indeed, an additional benefit of the use of cpn is represented by the possibility to realize multi scale models, but with the advantage of using the same modeling technique and software to combine, using hierarchical petri nets, intra cellular models and system models. 
continuous models, mainly based on differential equations, usually allow to qualitatively study the system but lack in description; conversely discrete models, such as agent based models and cellular automata, permit to describe in detail entities properties at the cost of losing most qualitative analyses. the methodology we present has advantages of both the two classical approaches based on continuous and discrete models, since it allows to gain good level of granularity in the description of cells behavior without losing the possibility of having a qualitative analysis. petri nets are a graphical modeling tool developed to model concurrency and synchronization in distributed systems. we propose a novel methodological approach that is based on high level pn, and in particular on colored petri nets, that can be used to model the immune system response at the cellular scale. computational modelingsystems biologypetri netsimmunityimmune systemdesign methodologyhumoral response15th international conference on bioinformatics incob 2016queenstown, singapore 21- september http://incob apbionet.org/issue-copyright-statement© the author 2016 furthermore, the presented methodology based on cpn allows the adoption of the same graphical modeling technique well known to life scientists that use pn for the modeling of signaling pathways. models and simulations allowed to study the immune system behavior, to test biological hypotheses about diseases and infection dynamics, and to improve and optimize novel and existing drugs and vaccines. mathematical and computational models showed to be a very important support tool for the comprehension of the immune system response against pathogens. to demonstrate the potentiality of the approach we provide a simple model of the humoral immune system response that is able of reproducing some of the most complex well-known features of the adaptive response like memory and specificity features. finally, such an approach may open the floodgates to the realization of multi scale models that integrate both signaling pathways models and cellular models built upon the same technique and software.
 mathematical and computational models showed to be a very important support tool for the comprehension of the immune system response against pathogens. models and simulations allowed to study the immune system behavior, to test biological hypotheses about diseases and infection dynamics, and to improve and optimize novel and existing drugs and vaccines. continuous models, mainly based on differential equations, usually allow to qualitatively study the system but lack in description; conversely discrete models, such as agent based models and cellular automata, permit to describe in detail entities properties at the cost of losing most qualitative analyses. petri nets are a graphical modeling tool developed to model concurrency and synchronization in distributed systems. their use has become increasingly marked also thanks to the introduction in the years of many features and extensions which lead to the born of “high level” pn.  we propose a novel methodological approach that is based on high level pn, and in particular on colored petri nets, that can be used to model the immune system response at the cellular scale. to demonstrate the potentiality of the approach we provide a simple model of the humoral immune system response that is able of reproducing some of the most complex well-known features of the adaptive response like memory and specificity features.  the methodology we present has advantages of both the two classical approaches based on continuous and discrete models, since it allows to gain good level of granularity in the description of cells behavior without losing the possibility of having a qualitative analysis. furthermore, the presented methodology based on cpn allows the adoption of the same graphical modeling technique well known to life scientists that use pn for the modeling of signaling pathways. finally, such an approach may open the floodgates to the realization of multi scale models that integrate both signaling pathways models and cellular models built upon the same technique and software. keywords computational modelingsystems biologypetri netsimmunityimmune systemdesign methodologyhumoral response15th international conference on bioinformatics incob 2016queenstown, singapore 21- september http://incob apbionet.org/issue-copyright-statement© the author 2016 