Blog - connections: Petri nets and beyond (part 1) (Rev #1)

I like knowing how things are connected. I hate when I find a new concept, sense that concept is connected to something I already know, but have no idea how to confirm my intuition. The people and writings at the Azimuth Project have introduced me to multiple new things. In this blog I explore one of these, Petri nets, and how they can build connections to some other islands in my knowledge.

Unlike many people here, I’m not good at math. Instead I depend on my software development skills. My approach here is to build a simple executable Petri net with software, and then convert it to other formats that can then be processed using tools from other domains. My Petri net example can be viewed as a mind map, processed by (at least) several Petri net tools, executed by various tools as a system of differential equations, displayed as a Unified Modeling Language (UML) class diagram and as a UML sequence diagram, converted into System Biology Markup Language (SBML) so it can be processed by numerous other tools, and run in a rectangular grid as an agent based model (ABM).

Petri nets have been well described elsewhere at the Azimuth Project, such as here, here, and here. Basically, a Petri net consists of a set of nodes and a set of edges that connect nodes. The nodes are of two types, **place** nodes (sometimes called states) and **transition** nodes. An edge is called an **arc**. Each place node contains zero or more tokens. Places only connect to transitions, and transitions only connect to places. Transitions act on the tokens in the places they connect to, by moving tokens from input places to output places.

Chemical reaction networks (CRN) were independently developed, but are almost exactly the same as Petri nets. They use different terminology and are typically used to represent chemical reactions. An excellent source of information on CRNs are Martin Feinberg’s lectures. The following diagram shows Feinberg’s first example CRN.

Xholon is a general-purpose Java-based modeling and simulation tool I’ve written. After being inspired by glowing accounts of Petri nets and CRNs at the Azimuth Project, I added a Petri net and CRN plug-in to Xholon. Xholon models are typically specified as XML hierarchical structures, with additional links between nodes.

The following is an example of a Petri net specified in Xholon XML format. It’s an implementation of the first example in Feinberg’s lectures. There are five places (A B C D E) each with a specified initial number of tokens. There are six transitions (A_BB BB_A AC_D D_AC D_BE BE_AC), each with a reaction rate represented as a symbol (a Greek letter) and a value k. Each transition has one or more input arcs and output arcs. Each arc has a weight and connects to one place node. The connectors specify the place as an XPath expression, which provides directions on how to get from the InputArc node to the place node. For example, *ancestor::PetriNet/Places/A* means locate an ancestor node called *PetriNet*, locate a child of that node called *Places*, and then a child of that node called *A*.

```
<?xml version="1.0" encoding="UTF-8"?>
<ReactionNetworkSystem>
<PetriNet roleName="Feinberg 1.1" kineticsType="1" p="1.0">
<QueueTransitions/>
<AnalysisPetriNet/>
<AnalysisCRN/>
<AnalysisCat/>
<Places>
<A token="140"/>
<B token="180"/>
<C token="200"/>
<D token="25"/>
<E token="80"/>
</Places>
<Transitions>
<A_BB k="0.0039" symbol="α">
<InputArcs>
<InputArc weight="1" connector="ancestor::PetriNet/Places/A"/>
</InputArcs>
<OutputArcs>
<OutputArc weight="2" connector="ancestor::PetriNet/Places/B"/>
</OutputArcs>
</A_BB>
<BB_A k="0.0038" symbol="β">
<InputArcs>
<InputArc weight="2" connector="ancestor::PetriNet/Places/B"/>
</InputArcs>
<OutputArcs>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/A"/>
</OutputArcs>
</BB_A>
<AC_D k="0.0037" symbol="γ">
<InputArcs>
<InputArc weight="1" connector="ancestor::PetriNet/Places/A"/>
<InputArc weight="1" connector="ancestor::PetriNet/Places/C"/>
</InputArcs>
<OutputArcs>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/D"/>
</OutputArcs>
</AC_D>
<D_AC k="0.0036" symbol="δ">
<InputArcs>
<InputArc weight="1" connector="ancestor::PetriNet/Places/D"/>
</InputArcs>
<OutputArcs>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/A"/>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/C"/>
</OutputArcs>
</D_AC>
<D_BE k="0.0035" symbol="ε">
<InputArcs>
<InputArc weight="1" connector="ancestor::PetriNet/Places/D"/>
</InputArcs>
<OutputArcs>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/B"/>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/E"/>
</OutputArcs>
</D_BE>
<BE_AC k="0.0034" symbol="ξ">
<InputArcs>
<InputArc weight="1" connector="ancestor::PetriNet/Places/B"/>
<InputArc weight="1" connector="ancestor::PetriNet/Places/E"/>
</InputArcs>
<OutputArcs>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/A"/>
<OutputArc weight="1" connector="ancestor::PetriNet/Places/C"/>
</OutputArcs>
</BE_AC>
</Transitions>
</PetriNet>
</ReactionNetworkSystem>
```

Xholon can transform this representation of the Petri net into Petri Net Markup Language (PNML). The PIPE2 Petri net tool displays it as the following, after some manual moving around of the nodes.

Diagram created by PIPE2 from Xholon data.

Xholon can also transform the model into Python syntax, so it can be executed by David Tanzer’s Petri net software.

```
# To run this Petri net (in Ubuntu Linux):
# (1) Save this python script to a .py file, for example pn_test.py .
# (2) Ensure that David Tanzer's petri1.py is in the same directory as this script.
# (3) Open a terminal window.
# (4) Run:
# python pn_test.py
# (5) Or run it while redirecting the output to a CSV file:
# python pn_test.py > pn_test.csv
from petri1 import *
# Petri net - Feinberg Lectures 1.1
A_BB = ("A_BB", [["A"]], [["B",2.0]])
BB_A = ("BB_A", [["B",2.0]], [["A"]])
AC_D = ("AC_D", [["A"], ["C"]], [["D"]])
D_AC = ("D_AC", [["D"]], [["A"], ["C"]])
D_BE = ("D_BE", [["D"]], [["B"], ["E"]])
BE_AC = ("BE_AC", [["B"], ["E"]], [["A"], ["C"]])
petriNet = PetriNet(
["A", "B", "C", "D", "E"], # states
[A_BB, BB_A, AC_D, D_AC, D_BE, BE_AC] # transitions
)
initialLabelling = {"A":140.0, "B":180.0, "C":200.0, "D":25.0, "E":80.0}
steps = 1200
petriNet.RunSimulation(steps, initialLabelling)
```

Petri nets turn out to be a really good starting point for modeling with other formalisms. For example, they can be automatically transformed into systems of differential equations. Let’s work through exactly what happens when we do this. If you don’t know or remember what differential equations are, then think of the following steps as a recipe, or driving instructions (how to get from Petrnet Town to Diffeq City), or maybe as a mathemagical incantation. This is something you can do by hand.

Write the name of each place (**A B C D E** in this model) on a separate line. Then write **d/dt ** in front of each name, and ** = ** after it. We say that **d/dt ** means *rate of change*, so **d/dt A** is the rate of change of the number of tokens in place A.

```
d/dt A =
d/dt B =
d/dt C =
d/dt D =
d/dt E =
```

For each place, look through all the transitions. When you find an input arc that involves that place, write ** -** at the end of the line. Write ** +** for an output arc. For example, **A** occurs five times in the six transtions, twice as an input and three times as an output. We say that **A** has five *terms*.

```
d/dt A = - + - + +
d/dt B = + - + -
d/dt C = - + +
d/dt D = + - -
d/dt E = + -
```

Each transition has a rate constant (symbolized as the Greek letters **α β γ δ ε ξ** in this model). For each term, write the rate constant directly after the **-** or **+** symbol. Each rate constant may also be given a numeric value (0.0039 0.0038 0.0037 0.0036 0.0035 0.0034 in this model), but in this process we will use the symbol. If you hate having to write Greek symbols, then write **alpha beta gamma delta epsilon xi** instead, for example **d/dt E = +epsilon -xi**.

```
d/dt A = -α +β -γ +δ +ξ
d/dt B = +α -β +ε -ξ
d/dt C = -γ +δ +ξ
d/dt D = +γ -δ -ε
d/dt E = +ε -ξ
```

Each transition has zero or more input arcs, where each arc references a place. For each term, write the names of all the places referenced by input arcs for that transition. If there is just one input arc, write the name of its place after the rate constant symbol. If there are two or more input arcs, write the names of all the places. These will all be multiplied with each other, so it’s OK to just write the single-character names of each place with no intervening spaces.

```
d/dt A = -αA +βBB -γAC +δD +ξBE
d/dt B = +αA -βBB +εD -ξBE
d/dt C = -γAC +δD +ξBE
d/dt D = +γAC -δD -εD
d/dt E = +εD -ξBE
```

If an input arc has a weight greater than 1, for example 2, then you can write it in one of two ways. For example, **BB** or **B^2**, because B times B is the same as B to the power of 2 (B squared). The symbol **^** is often used to mean *to the power of*.

```
d/dt A = -αA +βB^2 -γAC +δD +ξBE
d/dt B = +αA -βB^2 +εD -ξBE
d/dt C = -γAC +δD +ξBE
d/dt D = +γAC -δD -εD
d/dt E = +εD -ξBE
```

This is probably the step that’s the least intuitively obvious. I have problems with it. Any term that involves a transition where either the input arcs or output arcs have a weight greater than 1, may need a final adjustment. For example, if we’re figuring out the rate of change of B, then any term involving a transition that involves multiple instances of B, will need an adjustment. There are two such terms. Once I really get this, I’ll try to write it up in a more meaningful way. For now, the differential equation just works if we do things this way.

```
d/dt A = -αA +βB^2 -γAC +δD +ξBE
d/dt B = +2αA -2βB^2 +εD -ξBE
d/dt C = -γAC +δD +ξBE
d/dt D = +γAC -δD -εD
d/dt E = +εD -ξBE
```

So. This is our system of differential equations. The process we’ve just worked through has been captured as an algorithm in a Java program. Java is a popular programming language. To see the algorithm, look at the getDifferentialEquations() and getTerm() methods in the AnalysisCRN Java class. Warning: If you’re not a programmer, then this may be even scarier than thinking about differential equations for those of us not very math inclined.

What good are differential equations? A simple thing we can do is to format them so they look the way that mathematicians and others expect to see them when displayed on a web page (or in a book). Xholon automatically generated the following using a LaTeX format:

$\frac{d}{d t} A = -\alpha A +\beta B^2 -\gamma AC +\delta D +\xi BE$

$\frac{d}{d t} B = +2\alpha A -2\beta B^2 +\epsilon D -\xi BE$

$\frac{d}{d t} C = -\gamma AC +\delta D +\xi BE$

$\frac{d}{d t} D = +\gamma AC -\delta D -\epsilon D$

$\frac{d}{d t} E = +\epsilon D -\xi BE$

There are many tools that know how to carry out the calculations specified by a system of differential equations. For example, if I write the six differential equations in this format, then I can give it to popular software programs such as *Matlab* and *GNU Octave* which will happily turn it into a nice chart showing how each place (A B C D E) changes over time. *GNU Octave* is one of many open-source (free) programs that are compatible with the very successful *Matlab* commercial product that’s used in laboratories and universities around the world. The .m file was produced by first transforming the Petri net into Systems Biology Markup Language (SBML) format, and then using a System Biology Format Converter project (SBFC) tool to convert it to Octave .m format. All of these transformations are done automatically by software.

Diagram created by QtOctave, a GUI for Octave, from SBFC-generated content.