Céline Kuttler

Cédric Lhoussaine

Mirabelle Nebut

Rule-based Modeling of Transcriptional Attenuation at the Tryptophan Operon

Transcriptional attenuation at E.coli's tryptophan operon is a prime example of RNA-mediated gene regulation. In this paper, we present a discrete stochastic model of the fine-grained control of attenuation, based on chemical reactions. Stochastic simulation of our model confirms results that were previously obtained by master or differential equations. Our approach is easier to understand than master equations, although mathematically well founded. It is compact due to rule schemas that define finite sets of chemical reactions. Object-centered languages based on the π-calculus would yield less intelligible models. Such languages are confined to binary interactions, whereas our model heavily relies on reaction rules with more than two reactants, in order to concisely capture the control of attenuation.

.

 analyze the sensitivity of attenuation through probability functions and, more generally, discrete master equations. This approach benefits from a rich probability theory that gives valuable insights and measurement capabilities. However, apart from rare exceptions, master equations can only be evaluated numerically, and not solved symbolically. Each biological system requires an ad-hoc master equation or probability function that is usually hard to design from the mechanistic intuition of the system.

 [12] allows direct stochastic simulation, starting from a given multiset of molecules and a set of chemical reactions. Rule-based models are intuitive in the sense that they describe molecular interactions and are simpler to modify and extend than models based on classical mathematical functions.

. However, rewriting n-ary to binary reactions is tedious and requires sufficient expressiveness of formal languages. Sequences of reactions need to be executed within atomic transactions, so that no other interactions intervene.

. We obtain a

Introduction

Transcriptional attenuation is a control mechanism deployed by bio-synthetic operons across bacterial species [START_REF] Gollnick | Complexity in regulation of tryptophan biosynthesis in Bacillus subtilis[END_REF][START_REF] Gutierrez-Preciado | New insights into regulation of the tryptophan biosynthetic operon in Gram-positive bacteria[END_REF][START_REF] Yanofsky | 1' mRNA15(t!1),Ribo(m!1),mRNA10(t), mRNA11(t),mRNA12(t),mRNA13(t),mRNA14(t),S2(s~bl) -> Ribo(m),S1(s~fr),S2(s~fr) @ 1.0 # rule schema 8 'RNAPto2' RNAP(d!1),DNA1(t!1),DNA2(t) -> RNAP(d!2),DNA1(t)[END_REF]. Operons are sequences of jointly transcribed genes, bio-synthetic operons encode enzymes for the synthesis of amino acids. Attenuation prematurely interrupts an ongoing round of the operon's transcription, in situations where the environment already contains a high concentration of the corresponding amino acid. Summarizing, it works as follows. First, the amino acid concentration determines the speed at which a ribosome translates the nascent messenger Rna (mRna). Second, the ribosome's position controls how the mRna folds into a two-dimensional structure. Finally, the mRna structure sets the end point of transcription.

Although attenuation has been investigated within bacterial systems since the 1970s [START_REF] Kasai | Regulation of the expression of the histidine operon in Salmonella typhimurium[END_REF][START_REF] Yanofsky | Attenuation in the control of expression of bacterial operons[END_REF], it attracted significantly less interest than the control of transcription initiation, that is mediated by Dna binding proteins. This changed in the 2000s after the discovery of regulatory mechanisms in higher organisms that exploit Rna properties [START_REF] Beisel | Design principles for riboswitch function[END_REF]. Quantitative investigations of Rna-mediated regulation gained momentum for therapeutic approaches and synthetic biology [START_REF] Barboric | A new paradigm in eukaryotic biology: HIV Tat and the control of transcriptional elongation[END_REF].

E.coli 's tryptophan (trp) operon is the best understood bio-synthetic operon. It allows the bacterium to synthesize the amino acid tryptophan upon need. Tryptophan regulation in E.coli relies on two further mechanisms beyond transcriptional attenuation, that are not considered in this paper. concise description by two ingredients, 13 rules schemas (introduced in Section 3) from which we generate our 71 chemical reactions, and n-ary chemical reactions.

By means of rule schemas, we represent finite sets of chemical reactions in a compact manner, which differ only in the choice of certain molecule parameters that are quantified over e.g. folding or binding state, or location. This idea is well known from logic programming [START_REF] Lloyd | Foundations of Logic Programming[END_REF], unification grammars [START_REF] Stuart | An Introduction to Unification-Based Approaches to Grammar[END_REF] or term rewriting [START_REF] Baader | Term rewriting and all that[END_REF].

N-ary reactions are indispensable to intelligible representations of the trp attenuator, as our work indicates. We hypothesize the same holds for many other biological cases. By n-ary reactions, we refer to rules with three or more inputs, as opposed to binary reactions. They allow to incorporate global control into models, of which our work distinguishes three categories.

First are conditions for rule application. Here, one among a rule's multiple inputs is neither consumed nor modified by the rule's application. However if this molecule wasn't available, the rule could not be applied. For instance, we use this mechanism to model that transcription only aborts if the nascent transcript has folded into the termination hairpin, i.e. the corresponding rule checks this later's presence.

In the second category, an actor undergoes a state change as a side effect of the reaction between two others. For instance, after translation has progressed beyond a certain threshold, a state change hinders the corresponding mRna components to form hairpins.

The third category allows to switch between abstraction levels: upon application of a rule, one actor is replaced by the enumeration of its individual constituents, or vice versa. We use this mechanism for dedicated control segments of mRna that can interact as a whole with other segments, or be processed stepwise. Depending on the circumstances our model opts for their representation either as a whole, or as the sequence of their constituents.

Paper outline. We review the biological background in Sect. 2, introduce our rule-based language in Sect. 3, and review related languages in Sect. 3.4. As a first example, we model the multi-step race between transcription and translation by four rule schemas in Sect. 4.1. We quantitatively investigate the hypersensitivity of this basic model by simulation with the Kappa Factory3 [START_REF] Vincent Danos | Scalable simulation of cellular signaling networks[END_REF]. The concurrent elongation model of Sect. 4.2 extends, it explicitly renders the simultaneity of transcription and translation of the same mRna. We use it to investigate the quantitative impact on the multi-step race when only part of the mRna is present at the beginning of the simulation. In Sect. 5, we present our model of the detailed attenuation mechanism at E.coli 's trp operon. We qualitatively reproduce and confirm results of Elf and Ehrenberg [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF] regarding the probability of uninterrupted transcription into the full operon as a function of the rate of trp-codon translation, and discuss the quantitative differences between our results and theirs.

Transcriptional Attenuation

Transcriptional attenuation prematurely interrupts a gene's ongoing transcription, or that of an operon, when the cell does not actually need the corresponding proteins. E.coli 's trp operon encodes enzymes for the biosynthesis of the amino acid tryptophan (Trp). Their production is attenuated if the bacterium's environment provides sufficient amounts of tryptophan to feed on.

In this section, we first review the principles of gene expression in bacteria. Then, we introduce ribosome-mediated transcriptional attenuation, a regulatory mechanism used across bacterial species, and how specifically it functions at E.coli 's tryptophan (trp) operon [START_REF] Gollnick | Trp operon and attenuation[END_REF][START_REF] Yanofsky | Attenuation in the control of expression of bacterial operons[END_REF][START_REF] Yanofsky | Transcription attenuation: once viewed as a novel regulatory strategy[END_REF]. We omit certain aspects documented in the biological literature that do not enter our formal model presented in this paper, such as the role of transfer Rna in translation, or the redundancy of the genetic code.

Transcription copies information content from a Dna sequence into an mRna sequence. It is carried out by an enzyme called RNA polymerase. Rnap initiates its work by binding to a distinguished short Dna sequence which indicates the beginning of a gene (or operon), from where Rnap starts assembling an mRna molecule. In the following elongation phase Rnap advances stepwise over Dna, extending the growing mRna nucleotide by nucleotide, progressing at an average rate of 50 nucleotides per second. Transcription terminates when Rnap encounters a terminator sequence on Dna, if an additional condition is then fulfilled. Figure 1 illustrates this additional condition, that depends on a property of the mRna being transcribed. The linear mRna sequence can fold into stable secondary structures, which due to their shape are called hairpins. In order for transcription to terminate, while the Rnap encounters a terminator sequence on Dna, the most recent portion of the transcript must be folded into a hairpin.

Translation reads out an mRna molecule into the corresponding sequence of amino acids. It initiates with the binding of a ribosome to the free end of an mRna, the other end of which is still being elongated by an Rnap. The ribosome advances along the mRna towards the Rnap in steps of codons, which are words of three mRna nucleotides. For each codon, the ribosome adds the corresponding amino acid to the growing sequence. These amino acid sequences later fold into three dimensional structures known as proteins. While the average rate of translation is 15 codons per second, each step of the ribosome is actually limited by the abundance of the currently required amino acid. The ribosome slows down on codons for which the corresponding amino acid is in short supply.

Transcriptional attenuation subtly couples the termination of an ongoing round of transcription to the translation efficiency of the first part of the nascent mRna, where this latter is limited by a critical amino acid (tryptophan, for the trp operon). The so-called leader sequence consists of the operon's first few dozen nucleotides. Attenuation boils down to a race between the Rnap transcribing the leader Dna, and the ribosome translating the leader mRna. In a nutshell, the attenuation race is as follows. If the amino acid of interest is abundant, the ribosome advances at its maximal speed, and the terminator hairpin forms. Transcription then aborts. Conversely if the critical amino acid is rare, the ribosome stalls early within the leader. The stalled ribosome inhibits the terminator hairpin, hence the Rnap wins the race, and transcription continues into the operon's protein coding regions.

Trp leader architecture. The leader's architecture is fundamental to attenuation at E.coli 's tryptophan operon. We distinguish four segments S 1 to S 4 within the leader mRna, see Fig. 2. Each pair of adjacent segments folds into a hairpin, if neither of the required segments is masked by a ribosome. Three different secondary structures can occur within the leader mRna of the trp operon. They are named by their respective roles in attenuation. The pairing of S 1 with S 2 represents the pause hairpin, that between S 2 and S 3 is the anti-terminator hairpin, and the terminator shown in Fig. 1 is the pairing of S 3 with S 4 .

Hairpin co-occurrence and mutual exclusion. Each leader segment can only participate in one hairpin at the same time. Most importantly, the anti-terminator prohibits the terminator hairpin by sequestering S 3 . Because both require S 2 , the pause hairpin excludes the anti-terminator. On the other hand, the pause hairpin (S 1 • S 2 ) and the terminator (S 3 • S 4 ) can co-occur, since they do not compete for a shared segment. With respect to our model, we need to mention that hairpin formation is fater than any other reaction in the system. It is also important to bear in mind that the segments become available one by one while Rnap transcribes the trp operon's leader, and that the leader transcript progressively forms hairpins whenever the ribosome's position allows. Hence, the leader mRna never indeed remains unfolded, as simplifyingly shown in Fig. 2.

The role of hairpins. The impact of hairpins on transcription is significant, they determine pausing of the Rnap and termination of transcription. As opposed to this, mRna hairpins do not impair translation. A translating ribosome disrupts hairpins along its way, without significantly slowing down. We now detail on the pause and terminator hairpins at E.coli 's trp operon.

Pause hairpin. After Rnap has transcribed the segments S 1 and S 2 , it remains stalled on a strong Dna pause site, while the mRna rapidly folds into the pause hairpin. This combination resembles the conditions for transcription termination, however, it is reversible: Rnap resumes transcriptions after a ribosome has arrived along the transcript and disrupted the pause hairpin. Let us consider the details of this initial configuration for the attenuation race in Fig. 3 (left). Rnap is stalled on the pause site on Dna, more precisely on the nucleotide that we refer to as Dna 0 . It has so far transcribed the leader up to and including its segments S 1 and S 2 , that have folded into the pause hairpin. The approaching ribosome disrupts the pause hairpin with its step onto the 7 th codon of the leader mRna, which is the first step from the initial conformation in Fig. 3. The attenuation race now starts.

Terminator hairpin. The Dna leader of the trp operon contains a terminator sequence, just after the portion that encodes S 4 , see Fig. 1. When Rnap arrives here, the terminator mRna hairpin can form. The combination of terminator mRna hairpin and terminator Dna sequence aborts transcription. However if the anti-terminator is already present when S 4 is completed -it can appear as soon as S 2 and S 3 have been transcribed -the terminator is prevented. In this case Rnap continues unhindered through the terminator sequence, and reaches the enzyme-coding region of the operon (the structural genes).

Trp codons within the leader mRNA. We have not yet mentioned the codons 10 and 11 of the leader mRna (see Fig. 2), the translation of which each requires one tryptophan molecule. These two control codons determine the outcome of attenuation race. They act as sensors for the tryptophan concentration, and determine the speed of the ribosome's forward movement.

If tryptophan is in rare, the ribosome stalls on the control codons, hence its footprint does not advance far enough to mask the second segment. Soon later, the anti-terminator hairpin forms between S 2 and S 3 , and transcription continues into the structural genes. This is depicted as read-through configuration in Fig. 3.

Conversely if tryptophan is abundant, the ribosome efficiently translates through the control codons. From the time point the ribosome has reached the Fig. 3. Starting point and possible outcomes of the attenuation race at E.coli 's trp operon. Initial conformation (left): Rnap is paused by the pause hairpin, awaiting to be released by the ribosome's next step. Readthrough conformation: when tryptophan supply is low, the ribosome stalls on the control codons, the anti-terminator hairpin forms and transcription continues into the operon. Termination conformation: when tryptophan supply is high, the ribosome rapidly translates over the control codons. Before it unbinds from the mRna, the terminator hairpin forms and transcription aborts. 13 th codon, and until it dissociates from the stop codon 15, the ribosome's footprint masks S 2 , which prevents the anti-terminator [START_REF] Roesser | Ribosome release modulates basal level expression of the trp operon of Escherichia coli[END_REF]. The ribosome's unbinding delay from the stop codon is generally one second, which is a considerably long time scale, compared to all other reactions in the system. While S 2 remains blocked by the ribosome, Rnap continues transcription, it completes S 3 when reaching the 36 th Dna nucleotide, and S 4 at the 47 th Dna nucleotide. The terminator hairpin then forms and transcription aborts -this is the termination configuration.

Basal read-through due to premature ribosome release. A third possible outcome of the race is not covered by Fig. 3. When tryptophan supply is high, the ribosome occasionally dissociates from the stop codon sooner than expected. In that case S 3 can already have been transcribed, but S 4 not yet. Hence S 1 , S 2 and S 3 are available at the same time. With equal probability, either the pause hairpin or the anti-terminator forms, and in case of the latter, transcription continues. This basal read-through of the operon has been experimentally observed for 10 -15 % of initiated transcripts when tryptophan is abundant [START_REF] Nakamura | Regulation of basal level expression of the tryptophan operon of Escherichia coli[END_REF].

Rule Schemas for Chemical Reactions

In this section we first provide a formal and minimal rule-based language tailored to our needs (Sect. 3.1). We define chemical reactions, that operate on multisets of complex molecules with attributes such as Rnap • Dna [START_REF] Lloyd | Foundations of Logic Programming[END_REF]. Herein, the infix operator • indicates a complex between Rnap that is bound within a Dna sequence, more precisely at the position 23 stated by the attribute value of the Dna nucleotide. Other attributes of molecules could be the compartment of a molecule, or information on its states, for instance folding or binding state.

In Sect. 3.2, we present a language of rules schemas, that allows to define finite sets of chemical reactions in a compact manner. Rule schemas are like chemical reactions, except that attribute values are now extended to expressions with variables. All variables are universally quantified over finite sets, such that a rule schema defines a finite set of reactions. An example of a complex molecule is the term Rnap • Dna(x + 1) where x is a variable with values in {0, . . . , 50}. We introduce our language's stochastic semantics in Sect. 3.3.

As discussed in Sect. 3.4, more general ruled-based languages might have been used for our modeling study. The language in this paper is not intended as a contribution on its own, for the sake of its simplicity we however chose it for our presentation. Indeed, we relied on the software tool for another rule-based language to implement our models of Sect. 4 and 5.

Chemical Reactions

In order to define the syntax of attributed molecules, we fix a possibly infinite set of attribute values C and a finite set N of molecule names. We assume that each molecule name N ∈ N has a fixed arity ar(N ) ≥ 0, which specifies its number of attributes.

Molecules M ∈ Mol ::= N (c1, . . . , cn) | M1 • M2 Solutions S ∈ Sol ::= M | S1 , S2 Reactions S1 → k S2 Table 1. Chemical reactions where N ∈ N , c 1 , . . . , c n ∈ C, ar(N ) = n and k ∈ R + ∪ {∞}. Expressions e ∈ Exp ::= x | c | f (e1, . . . , en) Schematic molecules M ∈ SMol ::= N (e1, . . . , en) | M1 • M2 Schematic solution S ∈ SSol ::= M | S1 , S2 Rule schema ∀x1 ∈ D1 . . . ∀xn ∈ Dn. S1 → k S2 where V(S1) ∪ V(S2) ⊆ {x1, . . . , xn} Table 2. Rule schemas where x ∈ V, c ∈ C, f ∈ F, e 1 , . . . , e n ∈ Exp, N ∈ N , ar(N ) = n, D 1 , . . ., D n ⊆ C are finite sets, and k ∈ R + ∪ {∞}.
A molecule M , defined in Table 1, is a complex of attributed molecules. We write M 1 • M 2 for the complex of M 1 and M 2 . For instance, if Rnap, Dna ∈ N and 47 ∈ C then Rnap • Dna(47) is a molecule complex consisting of an Rnap that is bound to the Dna nucleotide at position 47. A chemical solution S is a multiset of molecules.

A chemical reaction is a rule that rewrites a solution S 1 into a solution S 2 , it is assigned a possibly infinite stochastic rate constant k ∈ R + ∪ {∞}. For instance, the following reaction states that an Rnap bound to the Dna nucleotide at position 23 may advance to the Dna nucleotide at position 24. The speed of this reaction is 50 sec -1 :

Rnap • Dna(23) , Dna(24) → 50 Dna(23) , Rnap • Dna(24)
In order to represent transcription, one would need many similar rules for the many other Dna nucleotides with different positions. This motivates the introduction of rule schemas, that allow to define such sets of chemical reactions in a compact manner.

Rule Schemas

In order to define rule schemas for chemical reactions, we introduce variables x for attribute values and expressions such as x + 1, in order to compute corresponding attribute values. By universal quantification over a finite set, we generalize the above chemical reaction to the following rule schema:

∀x ∈ {1, . . . , 49}. Rnap • Dna(x) , Dna(x + 1) → 50 Dna(x) , Rnap • Dna(x + 1) (0)
We thus need a set V of variables that are ranged over by x, and a finite set F of function symbols f ∈ F with arities ar(f ) ≥ 0. Furthermore, we assume an interpretation f : C ar(f ) → C for every f ∈ F. An expression e with values in C is a term with the abstract syntax given in Table 2. In our modeling case studies, we will assume that symbol + ∈ F of arity two is interpreted as addition on natural numbers. We freely use infix syntax as usual, i.e. we write e 1 + e 2 instead of +(e 1 , e 2 ). Given a variable assignment α : V → C, every expression e ∈ Exp denotes an element e α ∈ C that we define as follows:

c α = c x α = α(x) f (e 1 , . . . , e n ) α = f ( e 1 α , . . . , e n α )
A schematic molecule M is like a molecule, except we now allow for expressions in attribute positions rather than attribute values only. A schematic solution S ∈ SSol is a multiset of schematic molecules. As usual, we write V(S) for the set of variables that occur in molecules of S. A rule schema specifies the domains of variables occurring in the schematic solutions of the rule by universal quantification over finite sets.

For every variable assignment α : V → C that maps variables to values in their domain, we can instantiate the rule schema to finitely many reactions. A schematic molecule M is mapped to a molecule M α ∈ Mol. Similarly, schematic solutions S ∈ SSol get instantiated to solutions S α ∈ Sol:

N (e 1 , . . . , e n ) α = N ( e 1 α , . . . , e n α ) M 1 • M 2 α = M 1 α • M 2 α S 1 , S 2 α = S 1 α , S 2 α
A rule schema is instantiated to a set of chemical reactions, by enumerating the chemical reactions for all variable assignments licensed by the quantifiers:

∀x 1 ∈ D 1 . . . ∀x n ∈ D n . S 1 → k S 2 = { S 1 α → k S 2 α | α : V → C, α(x 1 ) ∈ D 1 , . . . , α(x n ) ∈ D n }

Stochastic Semantics and Simulation

For the sake of completeness, we recall the stochastic semantics of chemical reactions and how to use them for the stochastic semantics with Gillespie's algorithm. This underlines that our biological modeling case studies are indeed expressed in a formal modeling language. The semantics of a set of chemical reactions is a continuous time Markov chain (ctmc). Note that, for modeling convenience, we allow infinite rate constant ∞. Chemical reactions with infinite rates always have the highest priority and are executed immediately, that is without time delay. Such extended ctmcs with infinite rate constants can actually be converted to regular ctmcs by elimination of immediate transitions, while preserving sojourn time (i.e. how long the Markov chain stays in a given state) and probability transitions (that is, given a current state, the probability to make a transition to another given state) 4 .

L ⊆ {1, . . . , n} ⊕i∈LMi ≡ S S → k S ′ ⊕ n i=1 Mi k -→ L S ′ , ⊕ i ∈L Mi r = P {(L,k)|S k -→ L S 1 ≡S ′ } k ¬∃L∃S ′′ .S ∞ -→ L S ′′ S r -→ S ′ n = ♯{L | S → S1 ≡ S ′ } m = ♯{L | S → S2} S ∞(n/m) -----→ S ′
Table 3. Stochastic semantics of chemical reactions with finite and infinite rate constants.

The states of the extended ctmcs are congruence classes [S] ≡ of chemical solutions S with respect to the least congruence relation ≡ that makes complexation and summation associative and commutative:

M 1 • M 2 ≡ M 2 • M 1 (M 1 • M 2 ) • M 3 ≡ M 1 • (M 2 • M 3 ) S 1 , S 2 ≡ S 2 , S 1 (S 1 , S 2 ) , S 3 ≡ S 1 , (S 2 , S 3 )
In ---→ S ′ where the corresponding probability is r = n/m. The number of occurrences of immediate reactions leading from S to a solution congruent to S ′ is n, and the number of all occurrences of immediate reactions starting from S is m.

Such transitions are invariant under structural congruence, i.e. for all

S 1 ≡ S ′ 1 and S 2 ≡ S ′ 2 it holds that S 1 r -→ S 2 (resp. S 1 ∞(r) ---→ S 2 ) if and only if S ′ 1 r -→ S ′ 2 (resp. S ′ 1 ∞(r) ---→ S ′ 2 ). We can thus define [S] ≡ r -→ [S ′ ] ≡ by S r -→ S ′ and [S] ≡ ∞(r) ---→ [S ′ ] ≡ by S ∞(r)
---→ S ′ as the transitions of the extended ctmc. Gillespie's algorithm for stochastic simulation takes as input a finite set of chemical reactions and a chemical solution S. If reactions with infinite rate constants are applicable, it computes n and m as defined above for each immediately reachable solution S ′ , and returns such an S ′ with probability n/m jointly with a null time delay. Otherwise, it computes the overall rate of all possible transitions R = {r|S r -→S1} r, returns with probability r/R a solution S 1 with transition S r -→ S 1 jointly with a time delay drawn randomly from the exponential distribution with rate r.

Language Design Choices and Related Rule-based Languages

Models with rule schemas are more compact than if only simple reactions were used, thus easier to read. Attributed molecules and expressions that manipulate them were introduced, in the context of biological modeling languages, in [START_REF] John | The attributed pi-calculus with priorities[END_REF]. Note that even if rule schemas could be defined solely by means of variables, function symbols allow a better control and precision of the collection of reactions that are generated. For example, without function symbols, we would need to resort to name sharing to represent Dna sequences 5 . Each Dna nucleotide would bear two parameters, one referring to its predecessor, the other to its successor. Given link names {ℓ 0 , . . . , ℓ 50 }, our previous rule (0) on page 9 reads as ∀x, y, z ∈ {ℓ 0 , . . . , ℓ 50 }.

Rnap • Dna(x, y) , Dna(y, z) → 50 Dna(x, y) , Rnap • Dna(y, z)

Then starting from a DNA sequence Dna(ℓ 0 , ℓ 1 ), Dna(ℓ 1 , ℓ 2 ), . . . , Dna(ℓ 49 , ℓ 50 ), this rule schema instantiates into more ground rules than needed. For example, the rule Rnap • Dna(ℓ 1 , ℓ 10 ) , Dna(ℓ 10 , ℓ 45 ) → 50 Dna(ℓ 1 , ℓ 10 ) , Rnap • Dna(ℓ 10 , ℓ 45 ) is a meaningless instance of the above schema. Indeed, it is never applicable if the above DNA sequence is not modified as it is expected in a correct model. All models written with our rule schemas can be compiled, by instantiation, to finite collections of simple and formally well-defined chemical reactions. Although reactions do not define a Turing-complete language [START_REF] Cardelli | On the computational power of biochemistry[END_REF], their expressiveness is sufficient for our purposes. Furthermore, such collections of reactions are supported by standard tools for stochastic simulation such as Dizzy [START_REF] Ramsey | Dizzy: stochastic simulation of large-scale genetic regulatory networks[END_REF] or the rule-based language BioCham [START_REF] Chabrier-Rivier | The biochemical abstract machine BioCham[END_REF].

Alternative rule-based languages with higher expressiveness are Turing complete, e.g. the graph rewriting language Kappa [START_REF] Vincent Danos | Rule-based modelling of cellular signalling[END_REF], BioNetGen [START_REF] Michael | Methods in Molecular Biology, chapter Rule-Based Modeling of Biochemical Systems with BioNetGen[END_REF], and bigraphs [START_REF] Krivine | Stochastic bigraphs[END_REF]. Their pattern based graph rewriting rules resemble schemas, but their semantics is not based on instantiations to ground rules. They rather directly apply to arbitrary subgraphs satisfying the pattern. In contrast to our approach, such patterns may describe infinitely many reactions. Furthermore, stochastic simulation is possible without inferring all those reactions on before hand. This generation process is uncritical in the present paper, since the overall number of reactions remains small, but is the bottleneck in other applications, where it grows exponentially [START_REF] Vincent Danos | Scalable simulation of cellular signaling networks[END_REF][START_REF] Yang | Kinetic monte carlo method for rule-based modeling of biochemical networks[END_REF]. Another promising language is LBS [START_REF] Pedersen | A language for biochemical systems[END_REF]. Its general purpose semantics allows for translations to different concrete semantics such as ODEs and ctmcs. LBS also features compact description of reactions with yet ∀i ∈{0,. . . ,n-1}.Rnap • Dnai , Dnai+1 →e 1 Dnai , Rnap • Dnai+1

(1)

Rnap • Dnan →∞ Rnap , Dnan (2)
∀i ∈{0,. . . ,m-1}.

Ribosome • mRnai , mRnai+1 →e 2 mRnai , Ribosome • mRnai+1 (3) Ribosome • mRnam →∞ Ribosome , mRnam (4)
Table 4. Rules for n transcription steps, in race with m translation steps.

another approach by means of parameterized modules, species expressions and "non-deterministic" species. These formal rule-based languages were designed and used so far to tackle protein-protein interactions that occur in cellular signaling such as metabolic pathways. In contrast to this, our rule-based model deals with a fine-grained mechanism of gene regulation.

Hyper-Sensitivity of Multi-Step Races

In this section we illustrate rule schemas for chemical reactions with a simple yet interesting example, borrowed from Elf and Ehrenberg (2005) [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF]. Abstracting away from its detailed control by mRna hairpins, transcriptional attenuation boils down to a plain race between the two competing multi-step processes of transcription and translation. As intuition easily confirms, the probability that transcription wins the race decreases as the ribosome speeds up, and vice versa. We present two rule-based models for this multi-step race.

The basic model of Sect. 4.1 investigates the hyper-sensitivity of attenuation depending on the respective number of transcription versus translation steps (n vs m). Using Elf and Ehrenberg's rate constants for transcription and translation, we reproduce the results of Fig. 2 in [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF].

In Sect. 4.2 we enrich our basic model by what we call concurrent elongation. An additional parameter m 0 denotes the number of codons contained in the initial solution, the remaining codons are dynamically spawn by the Rnap at simulation time. We show the impact of this additional level of concurrency, with respect to the attenuation race, through simulation.

Basic Model of Transcription and Translation

Elf and Ehrenberg demonstrated that the relative change in the probability that transcription wins the race can be much sharper, than the relative change in the ribosome's speed. As our work confirms, this hyper-sensitivity of attenuation is determined by the number of transcription steps (n) versus translation steps (m). We give a basic rule-based model that allows to reproduce the results of Elf and Ehrenberg. As we believe, our framework is easier to understand and less prone to error than the master equation approach, while compact and mathematically well-founded.

Model. The following initial solution describes the starting point of the multistep race, where the Rnap and ribosome are bound to the first positions of Dna and mRna respectively:

Rnap • Dna 0 , ⊕ n i=1 Dna i , Ribosome • mRna 0 , ⊕ m i=1 mRna i
We use the following notational conventions. Molecule names are N = {Rnap, mRna, Dna, Ribosome}, attribute values C = N 0 , function symbols F = {+}, variables V = {i}, and value parameters n, m ∈ C. Because our model's attributed molecules bear only few arguments, for the sake of presentation we slightly differ from the formal syntax introduced in Sect. 3. We write attributes as indices for molecule names, instead of parenthesizing them, e.g. Dna i instead of Dna(i). Moreover we write ⊕ n i=1 Dna i instead of Dna(1), . . . , Dna(n). Finally, we emphasize that each mRna i denotes one codon, which biologically speaking is a sequence of three individual mRna nucleotides, that the ribosome reads out in one step.

Our model's rule schemas are listed in Table 4. Rule 1 for the n steps of the transcribing Rnap from one Dna nucleotide to the next remains as in Sect. 3.2, where it was the running example. The translation rule 3 is analogous and reflects the ribosome's m steps over codons. The remaining rules 2 and 4 model the dissociation from the respectively last positions of Dna and mRna. Note that, bearing rate ∞, dissociation occurs without advance of the simulation clock. Hence it does not quantitatively affect our simulation results compared to the model of Elf and Ehrenberg, that does not include dissociation. In order to incorporate the control conditions at the trp operon, we will refine the dissociation rules in Sect. 5.

Simulation. The plot reporting our simulation results in Fig. 4 is organized as follows. The y-axis gives the probability that Rnap wins the race, on a scale between zero and one. It corresponds to the proportion of simulations in which Rnap dissociates before the ribosome does. The x-axis reports the translation rate on a logarithmic scale, that we vary from 0.01 to 100 codons per second in our simulations with the Kappa Factory [START_REF] Vincent Danos | Scalable simulation of cellular signaling networks[END_REF].

The three models that each contribute one curve in the plot only differ in their numbers of translation (m) versus transcription (n) steps. We combined (m = 1 vs n = 1), (m = 1 vs n = 50), and (m = 10 vs n = 50).

Let us compare the sensitivity of these three models. When m = 1 and n = 1, the probability curve decreases gently, already showing some non-linearity. Increasing the number of transcription steps to n = 50 steepens the curve, i.e. increases the sensitivity. The transition becomes even sharper when the number of translation steps reaches higher values (m = 10, n = 50). Such values hold for systems where, unlike at E.coli 's trp operon, attenuation is the sole control mechanism [START_REF] Kasai | Regulation of the expression of the histidine operon in Salmonella typhimurium[END_REF]. It is worthwhile pointing out that in this model each of the m translation steps potentially slows down the ribosome's advance, while at E.coli 's trp operon, only 2 in 9 steps do.

Our rate constants are calculated as in [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF], to ensure the outcomes of the three races are comparable. We keep the total time to perform the series of n transcription steps constant, such that 1/e ′ 1 = 1 sec. Thus, the rate constant for one individual transcription step out of n is e 1 = n • e ′ 1 . For one translation step the rate constant is e 2 = m • e ′ 2 . Hereby e ′ 2 is the average rate for m translation steps, which varies logarithmically between 0.01 and 100.

As the next section will show, our model can smoothly be extended by additional concurrent issues, that are more difficult to handle within the master equation approach.

Concurrent Elongation of mRNA

In the basic model, the multi-step race was represented by a ribosome and an Rnap advancing along two independent strands of mRna and Dna. Here we add what we call concurrent elongation to the multi-step race. The idea is to reflect that Rnap still elongates a transcript when the ribosome starts translating its older end. Translation can now become limited by the slower transcription: the ribosome can only translate those codons that have previously been produced by the Rnap. Our simulation results demonstrate that the outcome of the race depends on the length of the initially available mRna. , and so forth respecting the pattern nil, nil, codon. If only one codon is part of the initial solution (m 0 = 1) we obtain:

Ribosome • mRna 0 , Rnap • Dna nil 0 ,
Rnap • Dna nil n →∞ Rnap , Dna nil n (3) ∀i ∈ {0, . . . , m -1}. Ribosome • mRnai , mRnai+1 →e 2 mRnai , Ribosome • mRnai+1 (4) Ribosome • mRnam →∞ Ribosome , mRnam (5)
Table 5. Rules for concurrent elongation (n steps), where the transcribing Rnap adds one new codon to the solution every three Dna nucleotides.

For the sake of simplicity we do not show the Dna position corresponding to mRna m . When m 0 = m + 1, the initial solution reduces to that of Sect. 4.1: In addition to the constraint on Dna nucleotide alternation, we assume that the initial solution contains m 0 codons, that the rule set will lead to the dynamic supply of additional mm 0 codons, such that the final solution shall contain m+1 codons (allowing for m translation steps), and that 1 ≤ m 0 ≤ m < 1 3 n. The ribosome's target codon mRna m corresponds to the Dna position 3 • (mm 0 + 1) -1. Beyond this, we assume that Rnap eventually reaches its own target, the n th position of Dna, without injecting additional codons to the solution.

Ribosome • mRna 0 , ⊕ m i=1 mRna i , Rnap • Dna nil 0 , Dna nil 1 , Dna codon
Table 5 lists our rule schemas. The rule 1 for one step of the Rnap, in which no codon is produced, resembles rule 1 of Sect. 4.1. It applies when leaving nucleotides of the form Dna nil i , whether or not the step leaving the next nucleotide yields a codon. Hence the quantification over x ∈ {codon, nil} for Dna x i+1 . The complementary rule 2 injects a new codon into the solution when the Rnap leaves a nucleotide of the form Dna codon i . The new codon's index is calculated from the current Dna position i and the initially available number of codons by the arithmetic expression m 0 -1 + (i + 1)/3. By doing so, we ensure that Dna 2 yields mRna m0 , Dna 5 yields mRna m0+1 , etc. up to the ribosome's target mRna m .

The rule for the Rnap's dissociation from Dna (3) only marginally differs from that of the previous subsection (in that the nucleotide bears the second attribute nil), and the ribosome advance and release rules (4 and 5) remain just the same. Simulation. We simulated our concurrent elongation model within the Kappa factory with several combinations of m 0 , m and n. Figure 6 shows the outcome of the race distinguishing whether only one codon is initially present (m 0 = 1), or all (m 0 = m), for the same number of translation (m = 15) and transcription steps (n = 50). When all codons are contained in the initial solution (m 0 = 15), the simulation results reduce to those of the basic model, whereas for m 0 = 1, the simulation curve shifts to the right, meaning that the probability that transcription wins the race increases. Indeed for each translation step, the ribosome's advance is potentially limited by the polymerase, that needs to add a further codon to the mRna. Hence, even if translation is efficient, the polymerase wins more often than for m 0 = 15. In our simulations we observed a lesser shift for m 0 = 10, not included in the plot.

Our simulations underline that the outcome of the multi-step race is parameterized not only by the n transcription steps and the m translation steps, but also by the number m 0 of initially available codons. This last parameter only appears when the model integrates concurrent elongation. We can now summarize our analysis of the multi-step race parameterized by m, n and m 0 , in terms of the shape of the curve that represents the probability that the polymerase wins the race:

-As pointed out by Elf and Ehrenberg [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF], the ratio of m to n determines the curve's slope. They are the key parameters of the hyper-sensitivity of ribosome-dependent transcriptional attenuation. -Varying m 0 shifts the curve. The polymerase's chance to win increases with m 0 , when m and n remain fixed, because m 0 constrains the ribosome's advance along mRna. As we observed, the shift increases with the difference between m and m 0 .

Incorporating concurrent elongation into our model was facilitated by our rule-based approach with arithmetic. It would have been more difficult with probability functions. In a model that includes concurrent elongation, the positions of the ribosome and the Rnap are not independent. The advance of the former is limited by that of the latter. This point was not considered in [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF].

Modeling Transcriptional Attenuation

This section presents our rule-based model of ribosome-mediated transcriptional attenuation at E.coli 's tryptophan operon. It refines our basic model of Sect. 4.1 in several points. The messenger Rna's representation dynamically grows while we simulate Rnap's advance, similarly as in the concurrent elongation model. But whereas in Sect. 4.2, we only used individual codons as building blocks of the transcript, the attenuation model also features mRna segments as a whole. Explicit representations of S 1 , S 2 , S 3 , and S 4 allow us to smoothly cover the dynamics of secondary structure formation, and incorporate the regulatory impact of hairpins on transcription. We make one notable exception to our all-in-one representation of mRna segments. Regarding S 1 , we switch between two different abstraction levels depending on the context, either representing it as a whole, or enumerating its codon sequence (⊕ 14 i=10 mRna i ). After introducing our attenuation model, we present simulation results in Sect. 5.2, and then explain the quantitative differences between our results and those of Elf and Ehrenberg [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF] in Sect. 5.3.

Rule Schemas

Table 6 provides the rule schemas of our detailed attenuation model. The notational conventions are based on those of our basic model of Sect. 4.1. We use molecule names N = {Rnap, mRna, Dna, Ribosome, S}, where Dna nucleotides are unary, attribute values C = N 0 ∪ {fr, bl, hp}, function symbols F = {+}, and variables V = {i, n, m, t, x}. Molecules with two attributes S x i represent segments of the mRna leader. Their lower index i ∈ {1, 2, 3, 4} denotes the segment's number, and the upper index x the segment's state which is among:

free (fr): available for hairpin formation, -blocked (bl): masked by the ribosome's footprint, -hairpin (hp): complexed into a hairpin with a neighboring segment. For instance, S hp 1 • S hp 2 denotes the pause hairpin.

The initial solution for our simulations reflects the starting configuration for the attenuation race, depicted in Fig. 3 on page 7:

Rnap • Dna 0 , ⊕ 50 i=1 Dna i , Ribosome • mRna 6 , mRna 7 , mRna 8 , mRna 9 , S hp 1 • S hp 2 , mRna 15 (0)
Rnap has transcribed the leader up to and including S 1 and S 2 , that are paired into the pause hairpin, and is paused on the zero-th Dna nucleotide, that is followed by a sequence of 50. The ribosome has initiated translation and is located on the 6th codon of the transcript leader. We explicitly render the codons 6 to 9, which precede the segment S 1 , and the stop codon 15, that is located between the segments S 1 and S 2 . In contrast, we do not render the codons preceding 6, since they do not matter to the attenuation race, and for the same reason we will not provide rules for the initiation of transcription and translation.

Hairpin formation is covered by rule schema 1, be it for the pause hairpin, the anti-terminator or the terminator. Because hairpin formation occurs on a much faster time scale than any other reaction, we approximate it with an infinite rate constant.

Translation rules (schemas 2 to 7 in Table 6). Rule schema 2 covers the bulk of translation steps, that do not have side effects, nor depend on tryptophan availability or other side conditions. It bears the reaction rate constant e 2 = 15s -1 , i.e. the ribosome makes 15 steps over mRna per second, in average. Rule schema 3 deals with the ribosome's step over the tryptophan codons within the leader, i.e. the control codons 10 and 11, where the distinct elongation rate constant e 3 holds. We will vary e 3 within ]0, 15]s -1 in our simulations, while e 2 remains fixed. Starting from our initial solution (the above equation 0) the next important event is melting the pause loop S hp 1 • S hp 2 , as the ribosome steps from mRna 6 to mRna 7 . Two points are worthwhile noting in rule 4's right part. First, obviously since the pause loop is melt, S 2 's state becomes free -and one could similarly expect a stage change at S 1 . But second and more importantly, instead of switching S 1 's state, we pass from the abstraction of the segment as a whole, to the enumeration of the codons ⊕ 14 i=10 mRna i that make it up. The sequence enumeration remains part of the solution as long as the ribosome's footprint partially covers the first segment, i.e. until it dissociates from the stop codon. This implicitly sequesters S 1 from hairpin formation -which would instantaneously occur through schema 1 if both S 1 and S 2 were around and in their free state.

For the ribosome's step from mRna 12 to mRna 13 , we introduce two rules with distinct preconditions. The common result of both is to reflect that the second segment gets masked by the ribosome's footprint, i.e. both rules produce S bl 2 . When S 2 is initially free, rule 5 applies. Otherwise, S 2 is paired into the anti-terminator hairpin, and rule 6 handles its melting through the ribosome's advance. When the ribosome dissociates from the stop codon (rule 7), S 1 reassembles, and S 2 unblocks.

Transcription rules. Rules 8 to 13 in Table 6 represent transcription. Rule schema 8 represents one step of Rnap in the simplest possible fashion, that was already discussed in Sect. 4.1 with the rate constant e 1 = 50s -1 . It applies to all Dna positions from 0 to 50 with a few exceptions that we discuss in the order they are applied, starting from our initial solution.

Transcription resumes at position Dna 0 (rule schema 9) after the pause hairpin has been disrupted. This is witnessed by S 2 being either free or blocked. Note that it would have been simpler to check the absence of the pause hairpin, but such negative tests are neither supported by the language used in this paper, nor by most current rule-based frameworks and tools, a notable exception to a certain extent is offered by [START_REF] Michael | Methods in Molecular Biology, chapter Rule-Based Modeling of Biochemical Systems with BioNetGen[END_REF].

The remaining rules deal with the creation of new mRna segments, and (anti)termination of transcription. When Rnap steps over to Dna 36 , the Rna segment S 3 is injected into the solution (reaction 10). S 4 follows at Dna 47 (reaction 11). Transcription terminates on Dna 47 provided there is a terminator hairpin, see S hp 3 • S hp 4 in reaction 12. If conversely the anti-terminator is present transcription proceeds to Dna 48 (rule 13) and continues transcription into the operon.

Summary: global control by n-ary rules. Finally we summarize our use of n-ary rules, that separates into three categories. Such n-ary rules can not be rendered in an intelligible fashion within object-centric approaches limited to binary interactions, namely π-calculus based modeling languages.

In the first category, we check whether the current solution fulfills a certain prerequisite, e.g. contains a certain molecule, or a certain molecule in a specific state. The rules for abortion versus continuation of transcription (rules 12 and 13) depend on which hairpin is around, terminator or anti-terminator. Other prerequisites we check are unfortunately less intuitive: sometimes one would prefer to impose negative conditions on rule application, which however is neither supported by our language, nor most other current rule-based frameworks. For instance, rule 9 resumes transcription if the pause hairpin is absent, which is the case if S 2 's state is free or blocked.

The second category are reactions that actually occur between two molecules, but entail the modification of a third. Examples are the rules for blocking S 2 as the ribosome proceeds to mRna 13 : rule 5 blocks S 2 , rule 6 disrupts the antiterminator hairpin at the same time, such that the states of both S 2 and S 3 change. The third category is the abstraction level switching for mRna segments, that is assembling the first segment from the codons 10 to 14 in rule 7, versus splitting it in rule 4. different pathways lead to the anti-terminator, each of them corresponds to a distinct curve, and a third curve sums up.

Simulation

The curve ribosome stalling corresponds to anti-terminator formation while the ribosome remains stalled on the control codons. This predominates when trpcodon translation is slow, becomes rarer as the translation efficiency increases, and drops below 1% when e 3 ≥ 9s -1 .

In a second pathway, the anti-terminator hairpin forms after the ribosome has released from the stop codon 15. This represents the basal read-through level of the trp operon (see page 8). The corresponding curve ribosome release in our simulation plot starts from zero, steadily increases with the rate of trpcodon translation, reaches its maximal level around 4.5% at a rate of trp-codon translation of 4 s -1 , and remains stable henceforth.

Our results shown in Fig. 7 qualitatively confirm those of Elf and Ehrenberg [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF]. However, even if the curves have the same shape, the asymptotic decrease of ribosome stalling toward 0 is less sharp in our case. While in their work the curves for ribosome release and ribosome stalling cross at a rate of 6s -1 , in ours they already do so at 5s -1 . Moreover our experiments predict a rate of basal transcription of slightly under 5% when trp-codon translation is efficient, where Elf and Ehrenberg predict 8%. . Transitions between configurations of the trp leader that our model yields, depending on the relative positions of the ribosome versus Rnap. We list the first segment in its blocked state for the sake of presentation, when the solution indeed enumerates its sequence.

Discussion

We believe that the differences between our quantitative results and those of [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF] are due to the greater level of detail rendered by our model. Elf and Ehrenberg reduce attenuation to a race between the ribosome and the polymerase, but their model does not make hairpins explicit. Instead, they infer transcription probabilities from the relative positions of both the Rnap and the ribosome. As opposed to this, our model explicitly renders hairpins. Table 7 summarizes how the configuration of the system evolves, assuming simulations start from our initial solution (equation 0), and the ribosome has disrupted the pause hairpin. For each relative position of the ribosome on mRna (columns) and the Rnap on Dna (rows), the table states which segments have been transcribed so far, and which hairpins have formed. Arrows between the table's cells corresponds to possible transitions of our solution. We distinguish the following positions of interest for the Rnap on Dna:

• between DNA nucleotide 0 and 35: only segments S 1 and S 2 are contained in the initial solution, • between 36 and 46: Rnap has completed segment S 3 , • on position 47: Rnap has injected S 4 to the solution. The ribosome's positions of interest on the mRna leader are:

• between codon 7 and 12: the ribosome's footprint has not yet reached S 2 ;

• between codon 13 and 15: the ribosome's footprint masks S 2 until dissociation from the mRna. • off: the ribosome has dissociated from the transcript, making segment S 2 newly available for hairpin formation.

Elf and Ehrenberg distinguish two cases of anti-termination (reproduced in our simulations in Fig. 7): anti-termination during ribosome stalling and antitermination after ribosome release.

Stochastic trajectories leading to ribosome stalling descend our table's column [START_REF] Ciocchetta | Bio-PEPA: a framework for modelling and analysis of biological systems[END_REF][START_REF] Daniel | A general method for numerically simulating the stochastic time evolution of coupled chemical reactions[END_REF]. The ribosome then remains between the codons 7 and 12, most likely stalling on the control (trp) codons within the first segment. After the polymerase has transcribed segment S 3 , the anti-terminator hairpin S 2 • S 3 forms. This henceforth excludes the terminator hairpin, and transcription continues into the structural genes. We believe that our model and that of Elf and Ehrenberg show the same behavior for this first pathway.

We identified the following key difference between the models regarding the second pathway to anti-termination. Elf and Ehrenberg assume that the probability of anti-termination after ribosome release is half the probability of reaching the configuration (Rnap ∈ [START_REF] Yanofsky | Attenuation in the control of expression of bacterial operons[END_REF]46], ribosome OFF). The corresponding cell is highlighted in light gray in Table 7. It is important to note that we separate it into two sub-cells A and B. We refine Elf and Ehrenberg's assumption as follows: stochastic trajectories leading to anti-termination after ribosome release pass through the sub-cell B, but never through the sub-cell A.

Careful consideration of the table explains our refinement. The highlighted cell can be reached from either its top or left neighbor. Coming from the left neighbor (Rnap ∈ [START_REF] Yanofsky | Attenuation in the control of expression of bacterial operons[END_REF]46], ribosome ∈ [START_REF] Gollnick | Trp operon and attenuation[END_REF][START_REF] Gutierrez-Preciado | New insights into regulation of the tryptophan biosynthetic operon in Gram-positive bacteria[END_REF]), it is equally likely to reach sub-cell A that entails termination, as to reach B that entails anti-termination. However, coming from the top neighbor (Rnap ∈ [0, 35], ribosome OFF) makes a difference. Because the ribosome has unbound the mRna before segment S 3 was completed, the pause hairpin S hp 1 • S hp 2 is already part of the solution. Hence the anti-terminator does not appear once S 3 is completed. Thus, descending the column ribosome off, the system always reaches configuration A, and transcription always terminates.

Based on our detailed model, and contradicting Elf and Ehrenberg, we claim that the state (Rnap ∈ [START_REF] Yanofsky | Attenuation in the control of expression of bacterial operons[END_REF]46], ribosome off) does not lead with equal probabilities to anti-termination by ribosome release versus termination. The downward transitions into sub-cell A, that contains the pause hairpin S hp 1 • S hp 2 and hence excludes the anti-terminator, increase the termination probability6 . This line of reasoning agrees with experimental knowledge on the impact of early ribosome release [START_REF] Roesser | Ribosome release modulates basal level expression of the trp operon of Escherichia coli[END_REF], and may explain why our experiments predict anti-termination less often than Elf and Ehrenberg's.

Conclusion

We have shown that rule-based modeling provides concise and elegant models for the fine-grained mechanism of transcriptional attenuation, a problem left open by previous work on discrete event modeling of the tryptophan operon [START_REF] Simão | Qualitative modelling of regulated metabolic pathways: application to the tryptophan biosynthesis in E .coli[END_REF]. The core ingredients for our model are rule schemas and n-ary chemical reactions. The importance of n-ary reactions renders representations of this case of genetic regulation in object-centered languages such as the stochastic pi-calculus [START_REF] Kuttler | A stochastic pi calculus for concurrent objects[END_REF][START_REF] Phillips | Efficient, correct simulation of biological processes in the stochastic pi-calculus[END_REF][START_REF] Regev | Computational Systems Biology: A Calculus for Biomolecular Knowledge[END_REF] inappropriate, in practice. We used the Kappa factory for stochastic simulation [START_REF] Vincent Danos | Scalable simulation of cellular signaling networks[END_REF], which provided us convenient analysis tools.

In our model we identified positions of individual nucleotides and codons within Dna and Rna sequences by numbers, abstracted over by variables, and addressed successors by simple arithmetic. This technique has its limitations when polymers become more complex than simple lists. Alternatively, we could assign names to molecular domains, and memorise those names in attribute values of adjacent molecules, similarly to Kappa. As shown in Sect. 3.4, the cost for such an alternative is that meaningless instances of rule schemas may be generated.

In future work, we plan to compute the exact probability that the ribosome dissociates before the segment S 3 is formed. This corresponds to the additional pathway that is not rendered by the model proposed by [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF] [START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF]. This would formally prove what we conjectured to be the source of the quantitative difference between the two models. We can indeed compute such a probability because there are finitely many pathways and all pathways are terminating (i.e. the system always reaches a configuration in which no rule is applicable). In other words, one can exhaustively unfold the underlying ctmc and thus compute the probability associated to each pathway (i.e. each branch of the ctmc).

Fig. 1 .

 1 Fig. 1. Transcription terminates if the most recent portion of mRna is folded into a hairpin, when Rnap reaches a terminator Dna sequence.

Fig. 2 .

 2 Fig. 2. Leader region of the trp mRna. Adjacent pairs of the four segments S 1 to S 4 fold into alternative hairpins. The anti-terminator hairpin S 2 •S 3 promotes transcription of the operon, whereas the terminator S 3 • S 4 aborts it.

  Fig. 3. Starting point and possible outcomes of the attenuation race at E.coli 's trp operon. Initial conformation (left): Rnap is paused by the pause hairpin, awaiting to be released by the ribosome's next step. Readthrough conformation: when tryptophan supply is low, the ribosome stalls on the control codons, the anti-terminator hairpin forms and transcription continues into the operon. Termination conformation: when tryptophan supply is high, the ribosome rapidly translates over the control codons. Before it unbinds from the mRna, the terminator hairpin forms and transcription aborts.Figure reproduced with permission from Elf and Ehrenberg (2005) [11].

Fig. 4 .

 4 Fig. 4. Probability that transcription wins in the basic model, as a function of the average translation rate e ′ 2 , for different numbers of translation (m) versus transcription (n) steps.

Fig. 5 .

 5 Fig. 5. General initial solution for the concurrent elongation model, containing an mRna of length m 0 and a Dna of length n. The Dna is composed such that upon simulation, every three steps of the Rnap one new codon is spawn; the final solution contains m codons.

2 ,nil 3 ,, Dna codon 5 ,Figure 5

 2355 Figure5illustrates the general case. In addition to the constraint on Dna nucleotide alternation, we assume that the initial solution contains m 0 codons, that the rule set will lead to the dynamic supply of additional mm 0 codons, such that the final solution shall contain m+1 codons (allowing for m translation steps), and that 1 ≤ m 0 ≤ m <1 3 n. The ribosome's target codon mRna m corresponds to the Dna position 3 • (mm 0 + 1) -1. Beyond this, we assume that Rnap eventually reaches its own target, the n th position of Dna, without injecting additional codons to the solution.Table5lists our rule schemas. The rule 1 for one step of the Rnap, in which no codon is produced, resembles rule 1 of Sect. 4.1. It applies when leaving nucleotides of the form Dna nil i , whether or not the step leaving the next nucleotide yields a codon. Hence the quantification over x ∈ {codon, nil} for Dna x i+1 . The complementary rule 2 injects a new codon into the solution when the Rnap leaves a nucleotide of the form Dna codon

Fig. 6 .

 6 Fig. 6. Probability that transcription wins in the concurrent elongation model, as a function of the average translation rate e ′ 2 , for different numbers of initially available codons (m 0 = 15 vs m 0 = 1), but the same number of transcription steps (n = 50) and translation steps (m = 15).

Figure 7 Fig. 7 .

 77 Figure7plots the relative transcription frequency (y-axis) against the rate of trp-codon translation e 3 (x-axis). For each value of e 3 from 0 to 15 s -1 in steps of one, we performed 5000 Gillespie simulations of our model. Recall that two

Table 3

 3 ′ , as many times as they apply for some index set L, provided that no immediate reaction can occur,

	-S	∞(r)

, we introduce transitions S k -→ L S ′ stating that S can be reduced to S ′ by applying a chemical reaction with rate constant k ∈ R + ∪{∞} to the subset of molecules in S with positions in L. Positions are the indices in multisets such as M 1 , . . . , M n that we also write as ⊕ n i=1 M i . We next introduce two transitions -S r -→ S ′ , where r ∈ R + sums up all rate constants of chemical reactions reducing S to S

This preliminary implementation of the κ-calculus was available to us as beta-testers, the web-based tool Cellucidate (http://cellucidate.com) is its successor.

For such an elimination procedure, see[START_REF] Kuttler | A stochastic pi calculus for concurrent objects[END_REF] and references therein.

This is actually how we implemented our model of Sect. 4 in the Kappa Factory.

Another effect is due to the anti-terminator hairpin melting as the ribosome moves on to mRna13 that is included in our model, which lowers the ribosome stalling probability.

Acknowledgements. The Master's project with Valerio Passini at the Microsoft Research -University of Trento Center for Computational and Systems Biology sharpened our view of the system from a biological perspective and lead us to a rule-based approach. The previous Master's project of Gil Payet (co-supervised with Denys Duchier) had confronted us with the limitations of object-based approaches to the representation of the complex dependencies of transcriptional attenuation. We thank Maude Pupin, who was the first to point us at the intricate regulatory mechanisms at E.coli 's tryptophan operon, and Joachim Niehren for his valuable coaching. Plectix BioSystems kindly made the Kappa Factory available to us, and gave us useful support while we carried out the experiments reported in this paper. Finally, we thank the CNRS for a sabbatical to Cédric Lhoussaine, and the Agence Nationale de Recherche for funding this work through a Jeunes Chercheurs grant (ANR BioSpace, 2009-2011).

Table 6. Rule schemas for hairpin formation (1), translation (2-7), and transcription [START_REF] Vincent Danos | Scalable simulation of cellular signaling networks[END_REF][START_REF] Vincent Danos | Rule-based modelling of cellular signalling[END_REF][START_REF] Dematté | The beta workbench: A tool to study the dynamics of biological systems[END_REF][START_REF] Elf | What makes ribosome-mediated trascriptional attenuation sensitive to amino acid limitation?[END_REF][START_REF] Daniel | A general method for numerically simulating the stochastic time evolution of coupled chemical reactions[END_REF][START_REF] Gollnick | Trp operon and attenuation[END_REF].

A Kappa Rules for Section 5

In the following we list the Kappa code of our detailed attenuation model in Section 5, as we implemented it in Kappa Factory version 12.2.0. in order to run stochastic simulations. This encoding remains as close as possible to our rulebased model, notably it does not represent Dna or mRna as interconnected chains, hence does not use name sharing. Some comments on the syntax of Kappa seem appropriate. For instance consider the following reaction produced by rule schema 1:

This reaction is given the name 'LoopS2S3(anti-terminator)' and an infinite rate @ $INF. The reactants S2(s ~fr) and S3(s ~fr) have an attribute s both with values fr. The reaction produces the molecules S2(s ~hp!1) and S3(s ~hp! 1), where the attribute s has the value hp. The modifier !1 indicates that the molecules form a complex, which is linked by the edge 1.

# rule schema 1 'LoopS2S3 (anti-terminator)' S2(s~fr),S3(s~fr) -> S2(s~hp!1),S3(s~hp!1) @ $INF 'LoopS3S4 (terminator)' S4(s~fr),S3(s~fr) -> S4(s~hp!1),S3(s~hp!1) @ $INF 'LoopS1S2 (pause)' S1(s~fr),S2(s~fr) -> S1(s~hp!1),S2(s~hp!1) @ $INF # rule schema 2 'RiboTo8' Ribo(m!1),mRNA7(t!1),mRNA8(t) -> Ribo(m!2),mRNA7(t),mRNA8(t!2) @ 15.0 ... 'RiboTo15' Ribo(m!1),mRNA14(t!1),mRNA15(t) -> Ribo(m!2),mRNA14(t),mRNA15(t!2) @ 15.0 # rule sche0ma 3 'RiboTo11_Trp' Ribo(m!1),mRNA10(t!1),mRNA11(t) -> Ribo(m!2),mRNA10(t),mRNA11(t!2) @ 1.0 'RiboTo12_Trp' Ribo(m!1),mRNA11(t!1),mRNA12(t) -> Ribo(m!2),mRNA11(t),mRNA12(t!2) @ 1.0 # rule schema 4 'RiboTo7_MeltS1S2' Ribo(m!1),mRNA6(t!1),mRNA7(t), S1(s~hp!2),S2(s~hp!2) -> Ribo(m!2),mRNA6(t),mRNA7(t!2),mRNA10(t),mRNA15(t), mRNA11(t), S2(s~fr),mRNA14(t),mRNA13(t),mRNA12(t) @ 15.0 # rule schema 5 'blockS2_ribo@13' Ribo(m!1),mRNA13(t!1),S2(s~fr) -> Ribo(m!1),mRNA13(t!1),S2(s~bl) @ $INF # rule schema 6 'meltS2S3_ribo@13' Ribo(m!1),mRNA13(t!1),S2(s~hp!1),S3(s~hp!1) -> Ribo(m!1),mRNA13(t!1),S2(s~bl), S3(s~fr) @ $INF