The Azimuth Project
Connections: Petri nets, systems biology, and beyond

Connections: Petri nets, systems biology, and beyond

Summary: A brief but systematic exploration of various types of networks, and how they’re really all the same. Starting with a simple chemical reaction network and corresponding Petri net, I demonstrate how to transform these into systems biology networks, systems of differential equations, Unified Modeling Language (UML) diagrams, mind maps, agent based models, and more. Along the way I introduce third-party tools that are able to process each of the resulting formats.

My thinking is - if I understand X, and X can be transformed into Y, then I should be able to understand Y. And then if someone mentions Z and that it’s related to Y, then I have a better idea of where it fits.

Part of the purpose of all this is that I understand new things better if I can approach them one small step at a time. I now have a basic understanding of Petri nets. I can use that to learn about stochastic differential equations, etc.

This tutorial works through one simple example to demonstrate the commonality/continuity between a large number of different ways that people use to understand the structure and behavior of the world around us. These include chemical reaction networks, Petri nets, differential equations, agent-based modeling, mind maps, membrane computing, UML, SBML, SBGN, . The intended audience includes scientists, engineers, programmers, and other technically literate nonexperts. No math knowledge is required.

This is the start of a table similar to Table 1 in Baez and Stay “Rosetta Stone” paper. Each approach included on this page uses different terminology for essentially the same things. Ports connect active objects with collections of passive objects (inputs or outputs).

 PassiveActiveContainerPorts (in,out)
CRNspeciesreactionn.a.complex (reactant complex,product complex)
Petri netplacetransitionn.a.arcs (input arcs,output arcs)
Mem Compobjectrulemembrane?
ABM?agentgrid?
SBMLspeciesreactioncompartment(listOfReactants,listOfProducts)
Xholonpassive objectactive objectcontainerports
Symmetric Monoidal Category?morphism?object (tensor products of species)

Contents

Chemical reaction networks

In a chemical reaction network (CRN), reactions transform complexes of chemical species into other complexes of chemical species. In the diagram below, A, B, C, D and E are types of species, while α, β, γ, δ, ε and ξ represent types of reactions. Reaction α transforms one instance of species A into two instances of species B. Reaction β does the reverse. In reaction γ, one A and one C become one D. ξ converts a B and an E into an A and a C. For each reaction, there’s a complex of species at the start of the arrow, and a different complex of species at the end of the arrow. Complexes include A, 2B and B+E.

<!-- Created with SVG-edit - http://svg-edit.googlecode.com/ --> CRN Chemical reaction network A 2B α β A+C D γ δ B+E ε ξ

CRN summary (as reported by Xholon AnalysisCRN)

The summary states that a chemical reaction network (CRN) consists of three sets, a set of species S, a set of complexes C, and a set of reactions R. Each reaction has a corresponding rate symbol, so for example ε is the rate symbol for reaction D→B+E.

CRN = {S, C, R}
S = {A, B, C, D, E}
C = {D, A+C, A, 2B, B+E}
R = {A→2B, 2B→A, A+C→D, D→A+C, D→B+E, B+E→A+C}
R (rate symbols) [α, β, γ, δ, ε, ξ]

Reaction rates

Each reaction proceeds at a possibly different rate. Reaction α might be 3.7 times faster at converting A into pairs of B (2B) than β is at converting pairs of B back into A. In this case, reaction rate α has a value of 3.7 and reaction rate β a rate of 1.0. Each type of reaction has its own rate.

The word “reaction” can refer to the type of reaction, to instances of that type, or to the reaction rate. Often, the reaction rate has a different name than the type. The reaction type with reaction rate α might be called A–>2B or ABB_ or r1 or whatever.

Kinetics

Each species in a CRN has a concentration c that may vary over time. Initially a CRN might have a concentration of 10 instances of A, 20 of B, 13 of C, 32 of D, and 17 of E. If the CRN is allowed to run in a simulator, the concentrations will typically change until they reach some final valies.

Each reaction has a rule or function or process that it runs in the simulator, that tells the reaction exactly what is should do to transform instances of its input species into instances of its output species. A kinetics is an assignment of a rule to each reaction. Numerous types of kinetics are possible, including Petri net kinetics, and mass action kinetics.

A simulator typically proceeds in units of time called time steps, where one time step might for example represent one second or one nanosecond or one year or whatever. If the kinetics and the rate constant and the species concentrations are all defined, then we might say that at time step 23 reaction D→B+E transforms species D into species B and species E at a rate of 1.2 instances per second.

Petri nets

In a Petri net, transitions transform the markings in places into markings in other places. In the diagram below, A, B, C, D and E are places, while A_BB, BB_A, AC_D, C_AC, D_BE and BE_AC are transitions. Each place has a marking that shows how many tokens are currently at that place. For example, place A currently has a marking of 140 tokens. Each transition has zero or more input arcs and zero or more output arcs. Each arc has a weight that tells the Petri net kinetics how to proceed. As an example, each time step transition A_BB transforms 1 token from place A into 2 tokens in place B, so a marking of [140, 180] in time step 0 would become a marking of [139, 182] in time step 1.

The structure of this Petri net is exactly the same as the structure of the earlier CRN. The diagram conventions and terminology are different, but the concept is the same. A CRN can readily be transformed into a Petri net, and vice versa.

Diagram created by PIPE2 from Xholon data.

Petri net - Feinberg 1.1

Petri net summary (as reported by Xholon AnalysisPetriNet)

places [A, B, C, D, E]
transitions [A_BB:A → B + B, BB_A:B + B → A, AC_D:A + C → D, D_AC:D → A + C, D_BE:D → B + E, BE_AC:B + E → A + C]

Mass action kinetics

In mass action kinetics, the rate of change is proportional to the product of the concentrations of all species in the input complex, to use CRN terminology. For example, the rate of change for reaction B+E→A+C is proportional to the current concentration of species B times the current concentration of species E. This product is multiplied by the rate constant ξ to obtain the actual instantaneous rate.

Membrane computing

A membrane system has a membrane structure consisting of hierarchically nested membranes (a rooted tree). In this case the structure is just a single membrane which contains five objects (A, B, C, D and E) and six rules (α, β, γ, δ, ε and ξ).

<!-- Created with SVG-edit - http://svg-edit.googlecode.com/ --> Membrane computing Chemical reaction network β: {B,B} -------- {A} γ: {A,C} -------- {D} δ: {D} -------- {A,C} ε: {D} -------- {B,E} ξ: {B,E} -------- {A,C} membrane A B C D E α: {A} ------- {B,B}

Xholon

Everything in Xholon is a node in a single tree. A node can be a passive object, an active object, a container, or any combination of these three. I wrote Xholon partly to have a tool that works the way I think. I’ve written a Xholon mechanism that enables modeling and simulation of Petri nets and chemical reaction networks. Xholon can export these models in a wide variety of formats including all of those contained on this page.

The Xholon specification for the chemical reaction network shown at the start of this page, is:

<?xml version="1.0" encoding="UTF-8"?>
<ReactionNetworkSystem>
  <PetriNet roleName="Feinberg 1.1" kineticsType="1" p="1.0" k="0.0039" timeStepMultiplier="1">
    <QueueTransitions shouldAct="true"/>
    <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" count="10" 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" count="10" 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" count="10" 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" count="10" 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" count="10" 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" count="10" 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>
  <ConstantlyStirredPot rows="80" cols="120" gridCellColor="000020" shouldPlacesMove="true"/>
</ReactionNetworkSystem>

Xholon summary (as reported by Xholon AnalysisCRN)

The species concentrations are given as real numbers and as integers. The rate constants are values chosen arbitrarily for simulation of each reaction R.

S (concentrations) [140.0, 180.0, 200.0, 25.0, 80.0]
S (concentrations) [140, 180, 200, 25, 80]
R (names) [A_BB, BB_A, AC_D, D_AC, D_BE, BE_AC]
R (rate constants) [0.0039, 0.0038, 0.0037, 0.0036, 0.0035, 0.0034]

Systems of differential equations

A system of differential equations is just a statement of how each species changes over time. One of the nice things about chemical reaction networks with mass action kinetics, is that you can derive the differential equations just by inspecting the network diagram. The following were automatically generated:

ddtA=αA+βB 2γAC+δD+ξBE \frac{d}{d t} A = -\alpha A +\beta B^2 -\gamma AC +\delta D +\xi BE
ddtB=+2αA2βB 2+ϵDξBE \frac{d}{d t} B = +2\alpha A -2\beta B^2 +\epsilon D -\xi BE
ddtC=γAC+δD+ξBE \frac{d}{d t} C = -\gamma AC +\delta D +\xi BE
ddtD=+γACδDϵD \frac{d}{d t} D = +\gamma AC -\delta D -\epsilon D
ddtE=+ϵDξBE \frac{d}{d t} E = +\epsilon D -\xi BE

Systems Biology Markup Language

The Systems Biology Markup Language (SBML) is a popular machine-readable (XML) format for representing and exchanging models, especially computational models of biological processes. Models are specified in terms of species, reactions, compartments, kinetics laws, species concentrations, and units. Does this terminology sound familiar? At least 230 tools support SBML.

Diagram created by COPASI, a biochemical network simulator, from Xholon-generated SBML content.

COPASI SBML - Feinberg 1.1

Systems Biology Graphical Notation

The Systems Biology Graphical Notation (SBGN) is a “visual notation for network diagrams in biology”, “an effort to standardize the graphical notation used in maps of biological processes”. SBGN is loosely connected with SBML, and involves many of the same people and organizations. SBGN offers three types of diagram. The following is a SBGN Process Description diagram.

Diagram created by SBGN-ED from Xholon data.

SBGN - Feinberg 1.1

System Biology Format Converter

The System Biology Format Converter project (SBFC) converts SBML to various other formats, currently Octave, XPP, Scilab, CellML, BioPax, and SBGN images.

XPPAUT

XPPAUT “is a tool for solving differential equations, difference equations, delay equations, functional equations, boundary value problems, and stochastic equations”.

An SBFC converter generated an XPPAUT compatible file from the Xholon-generated SBML content. COPASI was also able to generate an XPPAUT compatible file.

GNU Octave

GNU Octave “is a high-level interpreted language, primarily intended for numerical computations. It provides capabilities for the numerical solution of linear and nonlinear problems, and for performing other numerical experiments. It also provides extensive graphics capabilities for data visualization and manipulation. … The Octave language is quite similar to Matlab so that most programs are easily portable”.

An SBFC converter generated an Octave compatible (.m) file from the Xholon-generated SBML content.

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

Octave - Feinberg 1.1

Mind map

Mind maps are a popular way to organize ideas and information hierarchically. Freemind is an open-source mind mapping tool, whose XML format is readable by various other mind mapping tools.

Diagram created by Freemind from Xholon data.

Freemind Mind Map - Feinberg 1.1

Unified Modeling Language - Class diagram

The Unified Modeling Language (UML) is TODO

Diagram created by yuml.me from Xholon data.

UML class diagram - Feinberg 1.1

Unified Modeling Language - Sequence diagram

TODO This isn’t a very good UML sequence diagram, but it does show that all six transitions are firing each time step, and that they fire in random order.

Diagram created by websequencediagrams.com from Xholon data.

UML sequence diagram - Feinberg 1.1

Agent Based Modeling

TODO The grid contains reaction/transition and species/place nodes that constantly move within the well-stirred pot. When reactions and species are co-located, they may interact with each other.

Agent Based Modeling - Feinberg 1.1

Master equation

TODO

I’m reading the Network Theory series, and finding part 6 to be especially interesting. I sort of get the idea, but I’d like to do an Agent Based Modeling (ABM) simulation of the various rabbit behaviors in part 6. The simulation would have 2 grids, a huge field with “an inexhaustible supply of rabbits … randomly roaming around”, and a cage for caught rabbits. A rabbit catcher would stand in “a certain area” of the field, would catch all rabbits during each interval of time, and place the rabbits in the cage. On average the catcher would catch one rabbit per time interval, sometimes one rabbit, sometimes no rabbits, and sometimes more than one rabbit. There would be one ABM behavior for each Petri net transition described in part 6: Catching rabbits, Dying rabbits, Breeding rabbits, Dueling rabbits, Kissing rabbits, and Brawling rabbits. Once I’ve done that, then I think I can come back to this page and write something appropriate, including some code that outputs a Petri net as a master equation.

Ideally I would want to write out a master equation in a format that could be loaded by some existing third-party tool that knows how to process them. This would confirm my understanding.

Feynman diagram

TODO

My thoughts on Feynman diagrams are similar to what I’ve written above about the master equation. They are similar to UML sequence diagrams, but I would like to be able to write Feynman diagrams in a format that can be read by an external Feynman diagram tool. If people familiar with Feynman diagrams agree that it looks valid, then it probably is one.

Category theory

TODO

Category summary (as reported by Xholon AnalysisCat)

objects [A, B, C, D, E]
morphisms [A_BB:A → B ⊗ B, BB_A:B ⊗ B → A, AC_D:A ⊗ C → D, D_AC:D → A ⊗ C, D_BE:D → B ⊗ E, BE_AC:B ⊗ E → A ⊗ C]

Category Theory has the concept of a functor (Rosetta Stone p.10), a mapping from one category to another. Do Petri nets and CRNs belong to the same category, or are they different categories, presumably subclasses of symmetric monoidal category? What about the other modeling approaches I mention in this tutorial? How does category theory relate to my pragmatic use of software to transform a model in one format into models in other formats that can then be loaded and executed by software written to support a different research area?

System dynamics

System dynamics is another approach to modeling. It has similarities with Petri nets, and uses differential equations. The 1972 Limits to Growth model used DYNAMO, the original system dynamics tool from Jay Forrester’s group. Today there are several popular tools, including STELLA and Vensim.

Xholon can generate a Vensim model, although the graphics don’t appear correctly because the file format is undocumented. The third-party SDM-Doc can document the resulting model.

TODO

  • time series line charts for various kinetics types