Design DNA
Design DNA
Design DNA
Supplementary data
"Data Supplement"
http://rsif.royalsocietypublishing.org/content/suppl/2012/01/03/rsif.2011.0800.DC1.htm
l
References
P<P
http://rsif.royalsocietypublishing.org/content/early/2012/01/03/rsif.2011.0800.full.html#
ref-list-1
Receive free email alerts when new articles cite this article - sign up in the box at the top
right-hand corner of the article or click here
Advance online articles have been peer reviewed and accepted for publication but have not yet appeared in
the paper journal (edited, typeset versions may be posted when available prior to final publication). Advance
online articles are citable and establish publication priority; they are indexed by PubMed from initial publication.
Citations to Advance online articles must include the digital object identifier (DOIs) and date of initial
publication.
J. R. Soc. Interface
doi:10.1098/rsif.2011.0800
Published online
Designing correct, robust DNA devices is difcult because of the many possibilities for
unwanted interference between molecules in the system. DNA strand displacement has been
proposed as a design paradigm for DNA devices, and the DNA strand displacement (DSD)
programming language has been developed as a means of formally programming and analysing these devices to check for unwanted interference. We demonstrate, for the rst time,
the use of probabilistic verication techniques to analyse the correctness, reliability and performance of DNA devices during the design phase. We use the probabilistic model checker
PRISM, in combination with the DSD language, to design and debug DNA strand
displacement components and to investigate their kinetics. We show how our techniques can
be used to identify design aws and to evaluate the merits of contrasting design decisions,
even on devices comprising relatively few inputs. We then demonstrate the use of these components to construct a DNA strand displacement device for approximate majority voting.
Finally, we discuss some of the challenges and possible directions for applying these methods
to more complex designs.
Keywords: DNA computing; formal verication; probabilistic model checking;
DNA strand displacement
1. INTRODUCTION
M. R. Lakin et al.
y
t*
x*
y*
t*
t
t*
x
x
x*
y*
t*
y
y
t*
x*
y*
t*
t*
x*
y
y
t*
x*
y*
y*
t*
y
t*
t*
x*
y*
t*
t*
x*
y*
t*
t*
x*
y*
t*
Figure 1. Toehold-mediated DNA branch migration and strand displacement. (Online version in colour.)
2. BACKGROUND
2.1. DNA strand displacement
DNA strand displacement [8] is a mechanism for performing computation with DNA molecules. Once initial
species of DNA are mixed together, strand displacement
systems proceed autonomously [9] as increases in entropy
(from releasing strands) and enthalpy (from forming
additional base pairs) drive the system forward [10].
These increases typically result from the conversion of
active gate structures into unreactive waste. Furthermore, because DNA strand displacement relies solely
on hybridization between complementary nucleotide
sequences to perform computational steps, these systems
require no additional enzymes or transcription machinery, which in turn allows experiments to be run using
simple laboratory equipment.
In most strand displacement schemes, populations of
single strands of DNA are interpreted as signals, whereas
double-stranded DNA complexes act as gates, mediating
changes in the signal populations. Within the system, the
computational mechanism is toehold-mediated branch
migration and strand displacement [11]. At the periphery
of the system, signal populations may be connected to
uorophores for human-readable output, or regulated
by custom-designed aptamer molecules that interface to
the biological environment. The latter example highlights a key strength of DNA-based computational
devices: the ability to interface directly with biological
systems [12,13].
Figure 1 presents example branch migration and
strand displacement reactions. Each letter in the gure
represents a distinct domain (a sequence of nucleotides)
and the asterisk operator (*) denotes the Watson
Crick (CG, TA) complement of a given domain.
J. R. Soc. Interface
Short domains (represented in colour) are known as toeholds, while long domains (in grey) are often referred to as
recognition domains. We assume that toeholds are sufciently short (410 nucleotides) that they hybridize
reversibly with their complements, whereas recognition
domains are sufciently long (greater than 20 nucleotides) to hybridize irreversibly [11]. Each single strand
is oriented from the 50 (left) end to the 30 (right) end,
and each double-stranded complex consists of hybridized
single strands with opposite orientations. We assume
that the underlying nucleotide sequences have been
chosen such that distinct domains do not interact at all.
In the rst reaction from gure 1, an incoming strand
binds to a gate because the t toehold domain in
the strand hybridizes with its exposed complement
in the gate, producing the intermediate complex on the
right-hand side. Because the incoming strand is only
held on by a toehold, this reaction can be reversed, causing the single strand to oat away into the solution. In the
second reaction, the x domain in the overhanging strand
matches the next domain along the double-stranded
backbone, which means that the branching point
between the overhanging strand and the backbone can
move back and forth in a random walk called a branch
migration. Eventually, the random walk may completely
detach the short x strand from the gate in a strand displacement. This reaction is considered irreversible
because the invading strand is now attached to the gate
by a long domain as well as a toehold. Note that if the recognition domain on the strand did not match the next
domain along the gate, then branch migration could
not proceed, and the incoming strand would eventually
unbind. We call such binding reactions unproductive.
The third reaction is another branch migration, though
in this case no strand is displaced because even after
the y domain has been displaced, the rightmost strand
is still attached by a toehold. The fourth reaction is a
(reversible) unbinding reaction in which the rightmost
strand spontaneously unbinds because of the low binding
strength of the toehold.
Binding, migration and unbinding reactions such as
those illustrated in gure 1 allow signal populations
M. R. Lakin et al. 3
http://lepton.research.microsoft.com/webdna
M. R. Lakin et al.
the expected number of phosphorylations before relocation occurs? For further details on probabilistic model
checking of CTMCs, see [19,20]. For a description of the
application of these techniques to the study of biological
systems, see Kwiatkowska et al. [21]. All of the models
discussed in this paper were analysed using PRISM [7], a
probabilistic model-checking tool developed at the
Universities of Birmingham and Oxford. Additional
details of probabilistic model checking using PRISM are
provided in 4.
3. RESULTS
In this section, we present a series of case studies
that demonstrate the application of probabilistic modelchecking techniques to the design of DNA strand
displacement systems. As mentioned previously, to do
so we have extended the Visual DSD tool [2] with the
capability to generate, from DSD designs, corresponding
model descriptions that can be directly analysed by the
PRISM probabilistic model checker [7]. We will study our
modied two-domain gate designs from 2.1, and present
analyses of various correctness, reliability and performance properties of the gates using PRISM. Finally, we will
construct an approximate majority voting system using
these components and show how additional approximation mechanisms can be adopted in order to verify
this system.
3.1. Transducer gates: correctness
We begin by considering one of the simplest reaction
gates: a transducer that takes an X signal as input
and produces a Y signal as output. The gate can be
thought of as implementing a unary chemical reaction
X ! Y. We will demonstrate the use of (non-probabilistic) model checking to debug strand displacement
systems, by detecting a bug in a awed design for a
transducer gate.
Our initial transducer design is specied by the
DSD code of gure 2. This also includes a denition
of single-stranded signals, where the S(N,x) module
denotes a population of N single-stranded signals carrying the X domain. We assume that the t^ toehold is
dened globally and shared by all gates and strands.
The T(N,x,y) module represents N parallel copies of
the transducer gate that implements the X ! Y reaction. The body of the module denition consists of
two populations of N complexes, and two populations
of N strands. The species present in the initial state
(when N 1) are shown in gure 3a. The rst gate
(b)
t
c.1
(1)
y
x
x*
(1)
t* y*
x
t* x*
(1)
c.1
(2)
c.1 t
t
(1)
c.1 t
t* c.1* a*
c.1* t* a*
c.1
t*
(1)
t* a*
(1)
M. R. Lakin et al. 5
(1)
(1)
c.1
t* x* t* c.1* a* t* a*
x
x*
t* y*
c.1 t
(1)
c.1* t* a* t*
(1)
Figure 3. (a) Initial species and (b) expected nal species for the transducer gate. (Online version in colour.)
by Visual DSD) that, when evaluated in a particular state of a model, return the number of reactive
strands and reactive gates in that state, respectively.
The variable output gives the number of output
strands (in this example, kt^ x2l) and N is the
number of parallel copies of the system. So, we say
that the execution was successful when all copies of
the gate have produced the required output and there
are no reactive gates or strands (other than output
strands) present.
Notice that, by denition, when the specication
all_done given above is true; no further reactions
can occur. Hence, such states of the model are deadlock
states (those with no outgoing transitions to other
states). We specify the correctness of the system
design by stating that: (i) any possible deadlock state
that can be reached must satisfy all_done and
(ii) there is at least one path through the system that
reaches a state satisfying all_done. These two
properties can be represented by formulae in the (nonprobabilistic) temporal logic CTL, which can be veried
by PRISM:
M. R. Lakin et al.
(a)
(b)
(1)
c.1
(1)
c.2
(1)
c.2
x2
(1)
x1
c.1 t
(1)
(1)
c.1
t* a*
c.1 t
x1
(1)
t*
x1
x2
c.2 t
c.2
(2)
c.2
t*
a*
t* x1* t* c.2* a*
x2
(1)
(1)
(1)
(1)
x0
(1)
x1
(1)
x0
c.1
t* x0* t* c.1* a*
t* a*
x1
c.2
t* x1* t* c.2* a*
x0
(2)
t* x0* t* c.1* a*
x0
c.1
x0
(2)
(1)
x0
x1
t* a*
c.1 t
t*
x1
x2
c.2 t
(1)
(1)
t*
(1)
(1)
Figure 4. Deadlock states for the two faulty transducers in series. (a) State 1 (error); (b) state 2 (success). (Online version in
colour.)
x0 t c.1 a t a
t* x0* t* c.1* a* t* a*
x0
c.1 a t a
t x0
t* x0* t* c.1* a* t* a*
c.1
t x2 c.2 t a
x1* t* x2* c.2*t* a* t*
x1
x2 c.2 t
a t
t x2 c.2
x1* t* x2* c.2*t* a* t*
x0
c.1 a t a
t x0
t* x0* t* c.1* a* t* a*
x1
a
t x0 t c.1 a
t* x0* t* c.1* a* t* a*
c.1
x1
a t
t x2 c.2
x1* t* x2* c.2*t* a* t*
x1
x2 c.2 t a t
x1* t* x2* c.2*t* a* t*
c.2
t x2
there are still some reactive species left at the end. The
ka t^l strand from the X0 ! X1 transducer prematurely
activates the X1 ! X2 transducer without producing the
Figure 5. Corrected transducer gate code, with an additional new a declaration that prevents crosstalk between
different gates.
probability
M. R. Lakin et al. 7
1.0
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
terminate
error
success
0.5
1.0
1.5 2.0
T (104 s)
2.5
3.0
Figure 6. Plot showing the probability for each possible outcome of the faulty transducer pair, after T seconds. (Online
version in colour.)
parameter T from the queries, we can use PRISM to compute the eventual probability of each outcome:
M. R. Lakin et al.
40
35
30
25
20
15
10
5
0
expected percentage
2
3
4
5
6
N (number of copies of transducer pair)
Figure 7. Plot showing expected percentage of leftover reactive gates in the nal state of the system against the
number of parallel buggy transducersthat is, the parameter
N in the system S(N,x0)| T(N,x0,x1)| T(N,x1,x2).
We observe that the expected number of unreactive (i.e. correctly executed) gates increases as we add more parallel
copies of the buggy transducer system. (Online version in
colour.)
Finally, in this section, we consider performance properties of DNA strand displacement systems, as opposed
to the reliability properties discussed earlier. Seelig &
Soloveichik [23] showed that the time for strand displacement circuits to execute increases linearly with the
depth of the circuit. We can verify such properties
computationally with our model of DNA strand displacement and PRISM. Using the corrected transducer
design from gure 5, we constructed DSD models of
the form S(1,x0)| T2(1,x0,x1)| . . .| T2(1,x
f k 2 1g,xk), for various values of k, which correspond
to transducer chains of increasing length. Note that we
only consider a single copy of every transducer in the
chain. We used the following PRISM temporal logic query
to compute the expected time to reach a state in which
all of the gates have nished executing:
4.5
4.0
3.5
3.0
2.5
2.0
1.5
1.0
0.5
0
2
3
4
5
6
K (size of transducer chain)
(b)
(1)
(2)
(1)
(1)
(1)
(1)
z
x
c
t
t
z
x* t* z*
(1)
(1)
c
(1)
c* t* y* t* a* t*
a t a
x t y t c
t* x* t* y* t* c* a* t* a*
M. R. Lakin et al. 9
(1)
x
t* x* t* y* t* c* a* t* a*
(1)
(1)
x* t* z* c* t* y* t* a* t*
(1)
Figure 10. (a) Initial species and (b) expected nal species for the catalyst gate C(1,x,y,z). Garbage collection results in only
inert structures being present among the nal species. (Online version in colour.)
BX !X X
b Y X ! X B
d B Y ! Y Y
10
(b)
(1)
(1)
z
t
c
z
t* z*
(1)
t
c
(2)
(1)
(1)
(1)
(1)
t
(1)
t x
t y
t c
t* x* t* y* t* c*
(1)
c* t* y* t* a* t*
x
t y
t c
t* x* t* y* t* c*
a t
(1)
a* t*
t* z*
a
(1)
a* t*
a
c* t* y* t* a* t*
(1)
Figure 11. (a) Initial species and (b) expected nal species for the alternative catalyst gate C_NoGC(1,x,y,z). Note that this
gate does not perform the additional garbage collection reactions that produced the completely inert structures seen in gure 10.
(Online version in colour.)
9.0
expected time (104 s)
8.0
7.0
6.0
5.0
4.0
3.0
2.0
1.0
0
2
4
6
N (number of copies of catalyst gate)
Figure 13. DNA strand displacement (DSD) code for a catalyst gate, which extends the C_NoGC gate from gure 9 by
using the constant keyword from the DSD language to
abstract away from population changes due to accumulation
of waste and depletion of fuel.
When an X and a Y meet, reactions (a) and (b) convert one of them into an auxiliary species B with equal
probability. Then, when a B meets an X or Y, reactions
(c) and (d ) convert the B to match the species that it
encountered. In this second step, the probability of a B
encountering an X as opposed to a Y depends on the
initial populations of X and Y in the system, and this
fact allows the system to amplify any excess population
of one species over the other to converge on a consensus
in which all of the species are converted either to X or to
Y. Furthermore, it was proved by Angluin et al. [25] that
the system converges with high probability to the population that was initially in the majority, if the original
margin is sufciently large.
Note that the above-mentioned reactions all involve
catalysts, like those from 3.3. Thus, we can use the catalyst gates discussed therein to implement this
system of chemical reactions in the DSD language, as
shown in the code in gure 13. In order to reduce the
number of reactions and to make model checking
more tractable, we will use a catalyst gate without
garbage collection. Even with this simplication, however, the fact that there are cycles in the chemical
reactions means that the system can potentially use all
of the available fuel, which causes the state space to
grow enormously. To counteract this effect, we modify
the gate designs from 3.3 further, using the constant
keyword of the DSD language. This keyword declares
that the population of a particular species should be
held constant across all reactions in the system, even
those reactions where it is produced or consumed. This
approximation can be used when the species is in
excess, allowing us to abstract away from depletion of
fuels and accumulation of waste, in cases when these
species are in very high concentrations. This helps us to
greatly restrict the size of the state space in the PRISM
model, by essentially collapsing any states that differ
only by their populations of fuels and/or waste products.
Because the population protocol is not guaranteed to
form a consensus around the species that was initially in
the majority, we use PRISM to compute the probabilities
of reaching each consensus state (X or Y ) in the DNA
J. R. Soc. Interface
probability of choosing X
M. R. Lakin et al.
11
1.0
0.9
1.0
0.8
0.7
0.5
4
5
6
7
8
9
10
0.6
0
5
0.5
4
init X
1 1
0.4
0.3
0.2
init X
Figure 14. Surface plot which shows the probability of reaching a consensus of X, for various initial populations of
X and Y. (Online version in colour.)
0.1
0
1.0
0.5
0.5
1.0
X1
X2
X3
X4
X5
Y1
Y2
Y3
Y4
Y5
0.5000
0.7468
0.8709
0.9341
0.9665
0.2531
0.5000
0.6843
0.8082
0.8869
0.1290
0.3156
0.5000
0.6537
0.7699
0.0658
0.1917
0.3462
0.5000
0.6349
0.0334
0.1131
0.2299
0.3651
0.5000
'
'
'
'
'
S1
L1
S1*
S1
L1
L1
R1
'
R1
12
N*
'
'
'
R
'
L1
*
N
S1*
S*
S*
N*
*
N
L
S
R
'
'
'
S2
S1
S2
S1*
S*
S2*
S1*
S*
S2*
R2
'
S
R2
L2
S1
S1*
S*
S1*
S*
S1
'
L2
'
'
S1
R2
L1
L1
R2
'
Figure 16. Elementary reduction rules of the DSD programming language. We let S denote the migration rate of a domain
sequence S and fst(S) denote its rst domain. We also let N+ and N 2 denote the binding and unbinding rates, respectively,
of a toehold N^. We assume that fst(R2) = fst(S2) for rule (RM). This ensures that branch migration is maximal along a
given sequence and that rules (RM) and (RD) are mutually exclusive. (Online version in colour.)
4. METHODS
4.1. Model simulation in DNA strand
displacement
The syntax of DSD described in 2 interprets systems as
well-mixed solutions: hence we impose a standard set
of structural equivalence rules (see Lakin & Phillips [2]
for more details). In addition, we assume that no
long domain and its complement are simultaneously
unbound anywhere in the system. This well-formedness
restriction ensures that two species can only interact
with each other via complementary toeholds.
The rules in gure 16 present reduction rules that
formalize DNA interactions in the DSD language. The
arrows are labelled with rates that are used to
J. R. Soc. Interface
13
t x
t* x* u*
t x
t* x*
u*
x
t* x* u*
M. R. Lakin et al.
t x u
t* x* u*
RB;N
RXk ;rk
t*
x*
y*
t*
t*
x*
y*
t*
PRISM
[7] is a probabilistic model checking tool developed at the Universities of Birmingham and Oxford.
It provides support for several types of probabilistic
models, including CTMCs, which we use here. Models
are specied in a simple, state-based language based
on guarded commands. Support for several other
high-level model description languages has also been
made available through language-level translations to
the PRISM modelling language. For example, PRISM has
the ability to import SBML [27] specications, which
have an underlying CTMC semantics. Translations
from stochastic process algebra such as PEPA and the
stochastic calculus [28] have also been developed.
Formally, letting R0 denote the set of non-negative
reals and AP denote a xed, nite set of atomic propositions used to label states with properties of interest,
a CTMC is a tuple (S,R,L) where:
PRISM
14
parameters
states
time (s)
N4
N5
N6
X 3,Y 5
X 4,Y 5
X 5,Y 5
57 188
284 641
1 160 292
240 286
674 066
1 785 250
4
26
145
266
860
2602
wide range of properties of individual circuit components constructing using this design. We showed
that the PRISM model checker can detect bugs owing
to crosstalk between gates, analyse quantitative properties such as reliability and performance, and compute
the probability of different possible outcomes of the
gates execution. These techniques can be particularly
useful during the initial stages of gate design. Even
model checking a single gate executing in isolation, as
in 3.2, can help us to identify errors in the design
that would be difcult to quantify using simulationbased methods. Although multiple simulation runs
can be used to approximate the probability of a given
error, performing large numbers of simulations can be
time-consuming, particularly in cases such as gure 7,
where the error probability is low for large number of
gates. More importantly, model checking can be used
to identify the source of an error, by providing a specic
execution trace of the behaviour that leads to its occurrence. As illustrated in 3.3, we can also use model
checking to investigate the potential of different
system designs, even when analysed using relatively
small numbers of inputs. Finally, we have shown how
applying additional abstractions to the populations of
fuel and waste species can allow us to scale up to verifying more complicated systems, such as the approximate
majority population protocol [25].
Nonetheless, model checking has its limitations. As
the species populations grow, the number of reaction
interleavings explodes, which causes problems for
naively scaling up to larger systems. Table 3 shows statistics for a selection of the largest models that we used
to generate the results in 3 (model checking was run
on a 2.80 GHz Dell PowerEdge R410 with 32 GB of
RAM). The table shows the size of each model and
the time required to check a single property. As
expected, model sizes grow rapidly as population sizes
are increased, meaning that models larger than those
shown in the table could not be analysed. In 3.4, we
had to approximate the populations of fuel and waste
species in the model as constant in order to prevent
the state space from becoming too large to generate.
This effect can be mitigated to an extent, for example,
by careful gate design: our modied two-domain gates
were deliberately constructed to minimize the number
of asynchronous steps required for garbage collection and sealing off used gates, which greatly expand
the state space. We also used a high level of abstraction (the Innite semantics of DSD [2]) to reduce the
number of reactions in the model as far as possible.
However, even with the cleverest gate design, the
M. R. Lakin et al.
15
REFERENCES
1 Phillips, A. & Cardelli, L. 2009 A programming language
for composable DNA circuits. J. R. Soc. Interface 6,
S419S436. (doi:10.1098/rsif.2009.0072.focus)
2 Lakin, M. R. & Phillips, A. 2011 Modelling, simulating
and verifying Turing-powerful strand displacement systems. In DNA computing and molecular programming:
17th Int. Conf., DNA 17, 19 23 September 2011 (eds
L. Cardelli & W. M. Shih), pp. 130 144. Pasadena, CA,
USA. Springer Lecture Notes in Computer Science, no.
6937. Berlin, Germany: Springer.
3 Duot, M., Kwiatkowska, M., Norman, G. & Parker, D.
2006 A formal analysis of Bluetooth device discovery.
Int. J. Softw. Tools Technol. Trans. 8, 621632. (doi:10.
1007/s10009-006-0014-x)
4 Steel, G. 2006 Formal analysis of pin block attacks.
Theor. Comput. Sci. 367, 257 270. (doi:10.1016/j.tcs.
2006.08.042)
5 Calder, M., Gilmore, S. & Hillston, J. 2006 Modelling the
inuence of RKIP on the ERK signalling pathway using
the stochastic process algebra PEPA. In Transactions on
16
10
11
12
13
14
15
16
17
18
19
20
21
J. R. Soc. Interface
22 Cardelli, L. In press. Two-domain DNA strand displacement. Math. Struct. Comput. Sci.
23 Seelig, G. & Soloveichik, D. 2009 Time-complexity of multilayered DNA strand displacement circuits. In DNA
computing and molecular programming: 15th Int. Conf.,
DNA 15, Fayetteville, AR, USA, 8 11 June 2009, Revised
Selected Papers, vol. 5877 (eds R. Deaton & A. Suyama).
Lecture Notes in Computer Science, pp. 144 153. Berlin,
Germany: Springer.
24 Soloveichik, D., Seelig, G. & Winfree, E. 2010 DNA as a universal substrate for chemical kinetics. Proc. Natl Acad. Sci.
USA 107, 53935398. (doi:10.1073/pnas.0909380107)
25 Angluin, D., Aspnes, J. & Eisenstat, D. 2008 A simple
population protocol for fast robust approximate majority.
Distrib. Comput. 21, 87 102. (doi:10.1007/s00446-0080059-z)
26 Zhang, D. Y. & Winfree, E. 2009 Control of DNA strand
displacement kinetics using toehold exchange. J. Am.
Chem. Soc. 131, 17 303 17 314.(doi:10.1021/ja906987s)
27 Hucka, M. et al. 2003 The systems biology markup
language (SBML): a medium for representation and
exchange of biochemical network models. Bioinformatics
9, 524531. (doi:10.1093/bioinformatics/btg015)
28 Priami, C., Regev, A., Shapiro, E. & Silverman, W. 2001
Application of a stochastic name-passing calculus to
representation and simulation of molecular processes. Inf.
Process. Lett. 80, 25 31. (doi:10.1016/S0020-0190
(01)00214-9)
29 Chandran, H., Gopalkrishnan, N., Phillips, A. & Reif,
J. H. Localized hybridization circuits. In DNA computing
and molecular programming, 17th Int. Conf. DNA 17,
Pasadena, CA, USA, 19 23 September 2011, vol. 6937
(eds L. Cardelli & W. M. Shih). Lecture Notes in Computer Science, pp. 6483. Berlin, Germany: Springer.
30 Sandilands, E., Akbarzadeh, S., Vecchione, A., McEwan,
D. G., Frame, M. C. & Heath, J. K. 2007 Src kinase modulates the activation, transport and signalling dynamics of
broblast growth factor receptors. EMBO Rep. 8, 1162
1169. (doi:10.1038/sj.embor.7401097)
31 Wolf, V., Goel, R., Mateescu, M. & Henzinger, T. 2010
Solving the chemical master equation using sliding
windows. BMC Syst. Biol. J. 4, 42. (doi:10.1186/17520509-4-42)
32 Henzinger, T., Mateescu, M., Mikeev, L. & Wolf, V.
2010 Hybrid numerical solution of the chemical master
equation. In Proc. 8th Int. Conf. on Computational
Methods in Systems Biology (CMSB10), pp. 5565.
New York, NY: ACM.
33 Goss, P. J. E. & Peccoud, J. 1998 Quantitative modeling of
stochastic systems in molecular biology by using stochastic
Petri nets. Proc. Natl Acad. Sci. USA 95, 67506755.
(doi:10.1073/pnas.95.12.6750)
34 Heiner, M., Gilbert, D. & Donaldson, R. 2008 Petri nets
for systems and synthetic biology. In Proc. SFM 2008,
vol. 5016 (eds M. Bernardo, P. Degano & G. Zavattaro),
Lecture Notes in Computer Science, pp. 215 264. Berlin,
Germany: Springer.
35 Christensen, S. & Petrucci, L. 2000 Modular analysis of
Petri nets. Comput. J. 43, 224 242. (doi:10.1093/
comjnl/43.3.224)