Modelling Approaches to Molecular Systems Biology

fear is the mind killer

List of Papers

This thesis is based on the following papers, which are referred to in the text by their Roman numerals.

I Fange, D., Nilsson, K., Tenson, T., Ehrenberg M. (2009) Drug efflux pump deficiency and drug target resistance masking in growing bacteria. PNAS, 106:8215-20

II Fange, D., Lovmar, M., Pavlov, M.Y., Ehrenberg, M. (2010) Not competitive enzyme inhibitors revisited. manuscript

III Fange, D. and Ehrenberg, M. (2010) Deleterious effects of fluctuations in parallel metabolic pathways. manuscript

IV Fange, D., Elf, J. (2006) Noise-induced Min phenotypes in E. coli.

PLoS Computational Biology 2:e80

V Fange, D., Berg, O.G., Sjöberg, P., Elf, J. (2010) Stochastic

reaction-diffusion kinetics in the microscopic limit. PNAS published ahead of print.

VI Hattne, J., Fange, D., Elf, J. (2005) Stochastic reaction-diffusion simulation with MesoRD. Bioinformatics 21:2923-4

VII Dennis, P.P., Ehrenberg, M., Fange, D., Bremer, H. (2009) Varying rate of RNA chain elongation during rrn transcription in Escherichia coli Journal of Bacteriology 191:3740-6

VIII Fange, D., Ehrenberg, M. (2010) Rate of ribosomal RNA transcription modulated by rrn operon sequence and RNA polymerase interaction.

manuscript

Reprints were made with permission from the publishers.

Contents

1 Introduction . . . . 9

2 Practical aspects of model analysis . . . . 11

2.1 Kinetics in spatially homogeneous systems . . . . 11

2.1.1 Mesoscopic description . . . . 12

2.1.2 The mesoscopic description in the macroscopic limit . . . . 13

2.2 Kinetics in spatially inhomogeneous systems . . . . 13

2.2.1 Mesoscopic description . . . . 14

2.2.2 The macroscopic limit of the mesoscopic description . . . . 15

2.2.3 The mesoscopic description in the microscopic limit . . . . 15

3 Antibiotics as enzyme inhibitors and growth rate modulators . . . . . 19

3.1 Drug target resistance masking . . . . 19

3.1.1 Antibiotic drug efficiency . . . . 19

3.1.2 Drug resistance . . . . 20

3.1.3 The model . . . . 20

3.1.4 Steady state growth rate bistability . . . . 22

3.1.5 Drug target resistance masking . . . . 23

3.1.6 Conclusions and outlook . . . . 24

3.2 In vitro characterisation of enzyme inhibitors . . . . 25

3.2.1 Rapidly equilibrating inhibitors . . . . 27

4 The effect of stochastic fluctuations on bacterial physiology and growth rate . . . . 29

4.1 Stochastic fluctuations in parallel anabolic pathways and their effects on cell fitness . . . . 29

4.1.1 The model . . . . 31

4.1.2 Stochastic fluctuations and the rate limiting pathways . . . . 31

4.1.3 Optimal growth rates . . . . 32

4.1.4 Conclusion and outlook . . . . 34

4.2 Noise affecting the phenotype of the Min-system . . . . 35

4.2.1 The model . . . . 36

4.2.2 The wildtype and the filamentous mutants . . . . 37

4.2.3 Noise induced phenotypes . . . . 37

4.2.4 Conclusions and outlook . . . . 39

4.3 Modulation of the transcription elongation rate on rrn operons in E. coli . . . . 40

4.3.1 The model . . . . 41

4.3.2 Initiation rates and operon transit times . . . . 41

4.3.3 Polymerase distribution . . . . 43 4.3.4 Conclusion and outlook . . . . 45 5 Summary in swedish – Systembiologisk modellering på molekylär

nivå . . . . 47 5.1 Dolda antibiotikaresistensmutationer . . . . 47 5.2 Stokastiska fluktuationer i aminosyrabiosyntes påverkar till-

växthastigheten . . . . 48 5.3 Stokastiska fluktuationer påverkar fenotypen hos

Min-proteinerna . . . . 48 5.4 Variation av transkriptionshastigheten på det ribosomala RNA

operonet . . . . 49

6 Acknowledgements . . . . 51

References . . . . 53

1. Introduction

Cells are complex creatures. While this is undoubtedly true for eukaryotic cells, which have organelles, clear internal structure, active transport, etc, etc, this is also true for the bacterial cell. The E. coli bacterial cell needs fine tuned spatial and temporal regulation of intracellular components. For example, the E. coli cell makes fine tuned decisions on when and where to divide, when to start DNA replication, and how many proteins should be made of each kind in order to utilise a nutrient efficiently.

We generally know a lot about bacteria such as E. coli

, especially con- cerning the different cell components, and how they interact. However, when it comes to questions of how the components interact together in space and time to form a phenotype the results may become harder to interpret. Here mathematical modelling can serve as a tool for rationalising the underlying mechanisms.

This thesis describes efforts to understand how the bacterial intracellular architecture gives rise to a certain phenotype, using mathematical models of the intracellular dynamics. The analysis is carried out for a number of different systems in E. coli, which are not connected in any other way than the need for mathematical models to understand a phenotype. Some effort have also been made to develop tools that make mathematical analysis possible.

Below follows an overview of the different parts of the thesis. Section 2 is a review of the different tools used in the model analysis in the different parts of the thesis. This section focuses on how to mathematically and numerically describe intracellular chemical kinetics. Its last part (2.2.3) describes the work in paper V, which makes the reaction-diffusion master equation possible to use also at high spatial resolution.

In section 3 the interactions between an antibiotic drug and its drug-target in the context of a living cell is described. In the first part (3.1) I explain the observed phenomenon of drug target resistance masking using a course grained model of bacterial growth in the presence of antibiotic drugs. Drug target resistance masking means that drug target affinity differences do not show up as growth differences. The model suggests a mechanism for drug target resistance masking, and also states the requirement for this to show up.

Antibiotic drugs are often enzyme inhibitors, and the second part (3.2) focuses on in vitro characterisation of enzyme-inhibitors. Here I give a short example

1Just the sight of the standard reference literature for E. coli [1], should make you certain of that

of how measurements of steady-state rate of product formation in the presence of enzyme inhibitors may lead to ambiguous interpretation of the data.

Section 4 describes three different systems where the discreteness of intra-

cellular molecules and the stochastic nature of chemical reactions affect the

phenotypes of the living cell. In the first part (4.1) I describe how stochastic

fluctuations in anabolic production pathways may affect the growth rate of a

bacterial cell. To do this, I have modelled amino acid turnover, and how this

is affected by stochastic fluctuations in the amino acid biosynthetic enzyme

pools. The second part (4.2) describes how stochastic fluctuations are essen-

tial for describing two mutant phenotypes of the Min system in E. coli. In a

mutant with changed membrane lipid composition the standard Min-protein

pole-to-pole oscillation pattern is changed into a pattern in which spots appear

and disappear on the membrane. The appearance of a spot is described as a

random nucleation event. Lastly, in the third part (4.3) I describe the results

of a model implemented to describe the RNA polymerase density variation on

the ribosomal RNA operon in E. coli, and also the reduced operon transit time

in repsonse to inactivation of 4 out of the 7 ribosomal operons.

2. Practical aspects of model analysis

To describe intracellular chemical kinetics one, in general, needs to consider both the discrete, stochastic changes due to chemical reactions and the spatial distribution of molecules in the cell [2]. However, when the spatial distribution is not crucial for the understanding of the system, a very common approxima- tion is to assume that the cell is well-stirred, and thus spatially homogeneous.

Another very common approximation is to use macroscopic rate equations.

In this approximation it is assumed that effects of stochastic fluctuations are negligible such that the number of molecules is well described by its average.

In this section I will go through the different levels of chemical kinetics that have been used in this thesis. I will start with the spatially homogeneous cases of chemical kinetics, and then move on to how chemical kinetics can be de- scribed in spatially non-homogeneous cases. In the latter section, I will intro- duce the result of paper V, which describes how to correct the reaction prob- abilities in spatially discrete systems such that they match the microscopic description (see section 2.2.3).

Before describing each method, I would like to point out that the best choice of method always depend on the problem at hand. The most complete descrip- tion must not always be the best description, as complex models are more time-consuming to work with, and usually harder to intuitively understand.

The best choice of method for the different problems described in this thesis is more (e.g. section 3.2), or less (e.g. section 3.1) obvious.

The sections below rely on examples of simple chemical reactions rather then on a complete, but abstract, description of chemical kinetics.

2.1 Kinetics in spatially homogeneous systems

In a well-stirred, spatially homogeneous system a complete description of its

state is given by the copy-number of each species [3, 4]. The well-stirred ap-

proximation implies that diffusion is fast enough to homogenise the spatial

distribution on the time scale of chemical reactions [4]. Chemical reactions are

here described as discrete, stochastic, memory-less, events changing the state

of the system. The time-evolution of the system can therefore be described by

a so called jump Markov process [4]. The loss of memory stems from the vast

number of non-reactive solute collisions that a reactant will make before the

next reaction event. This leads to randomisation of molecule-specific features, such as velocity and position [4].

Consider the production and consumption of X-molecules according to the following scheme

x

x + i x

x − j

(2.1)

Here g(x) and h(x) are the production and consumption rates of X-molecules, respectively, and i and j are the corresponding stoichiometries of the reactions.

When the system evolves as a jump Markov process, g(x) and h(x) depend only on the state of the system and not on time.

2.1.1 Mesoscopic description

Given that the system in eq. 2.1 is a jump Markov process, the probability, P(x,t + ∆t), of having x number of molecules at time t + ∆t is given by

P(x,t + ∆t) = P(x,t)

+ ∆t {g(x − i)P(x − i,t) + h(x + j)P(x + j,t)}

− ∆t {(g(x) + h(x)) P(x,t)}

(2.2)

where ∆tg(x − i) is the probability that the system changes from state x − i to x during the short time ∆t. Re-arranging and going in the limit ∆t → 0 results in the so called master equation [4, 3]

dP(x,t)

dt = g(x−i)P(x−i,t)+h(x+ j)P(x+ j,t)−(g(x) + h(x)) P(x,t) (2.3) The time-evolution of the X-molecule distribution is given by an infinitely large set of ordinary differential equations, since the number of X-molecules ranges from x = 0 to x = ∞.

The master equation can generally not be solved analytically [4], or even numerical solutions can be impossible. The reason for the latter is that even though the number of ODEs described by eq. 2.3 may be truncated at N and solved numerically, the case of a k-molecule-species requires N

number of ODEs. The other method is to approximate properties of the probability dis- tribution by making stochastic simulations.

Trajectories of the time-evolution of the system described by the master

equation in eq. 2.3 can be implemented using a Monte Carlo method known

as the stochastic simulation algorithm or the Gillespie algorithm [5]. The algo-

rithm is straightforward: randomly sample a reaction-time and a reaction from

their corresponding distributions. Update the system according to the sampled

time, and sampled reaction, and repeat the procedure. The probability of leav-

ing state x during the short time-interval ∆t is given by ∆t (g(x) + h(x)). In the limit ∆t → 0, the time-evolution of the probability of still having x number of X-molecules is given by

dP(x,t)

dt = − (g(x) + h(x)) P(x,t) (2.4)

and hence, the reaction-times are exponentially distributed with a mean time τ = (g(x) + h(x))

. The probability of the reaction associated with g(x) oc- curring is P

= g(x)/(g(x) + h(x).

Stochastic simulations of the master equation have been used for analysis of the systems in section 4.1 for the dynamics of mRNA and protein abundances, and in section 4.3 for the RNA polymerase translocation dynamics.

2.1.2 The mesoscopic description in the macroscopic limit

The mesoscopic description can be connected to a macroscopic description using the example below. By taking the average of eq. 2.3, or just by mass action, the average number of X-molecules, hxi , must follow [4, 6]

d hxi

dt = i hg(x)i − j hh(x)i (2.5)

When the numbers of X-molecules are very large, the stochastic fluctuations can be neglected, such that hxi ≈ x and hg(x)i ≈ g (hxi) = g(x) and thus

dx

dt = ig(x) − jh(x) (2.6)

This is the macroscopic description of chemical kinetics. The macroscopic description can be a very good approximation in vitro, where the number of molecules can be kept high by virtue of experimental design. It should, how- ever, be noted that the macroscopic equation in 2.6 does not even describe the average number of X-molecules, hxi, in an in vivo situation, where the number of molecules is low, and g(x) and r(x) can be non-linear, such that hg(x)i 6= g (hxi).

The spatially homogeneous macroscopic description has been used in sec- tion 3.1 to describe the intracellular dynamics of antibiotic drug and ribosome concentrations, and also in section 3.2 to describe the dynamics of in vitro experiments.

2.2 Kinetics in spatially inhomogeneous systems

When the conditions for spatial homogeneity described in the previous section

cannot be fulfilled, the positions of the molecules need to be considered. In

bacteria there are obvious cases in which the spatial aspect must be taken into account, such as the positioning of the septum in E. coli cell division, as studied in paper IV. There are also other, less obvious, cases such as the MAP- kinase phosphorylation cycle described in Takahashi et al. [7]. Here rebinding due to spatial proximity can change the overall behaviour of the system.

Chemical kinetics in spatially inhomogeneous systems can be described by two different methods. Either by a mesoscopic description (2.2.1) in which the physical space is discretised [6], or by a microscopic description where both space and time are continuous [8, 9].

2.2.1 Mesoscopic description

When diffusion is too slow to homogenise the spatial distribution of mole- cules in the cell in between reactive events, spatial homogeneity may still be achieved on shorter length scales. This lays the foundation for the reaction- diffusion master equation (RDME) [6]. In the RDME description, physical space is divided into subvolumes, and the state of the system is extended to x = {x

}, where x

is the number of X-molecules in subvolume n. In the RDME the state changes both by chemical reaction events in each subvolume, as in section 2.1.1 above, and by diffusion events, where one molecule jumps from subvolume n to a neighbouring subvolume m. Consider again the reaction scheme in eq. 2.1. In the RDME description, the probability P(x

,t) of having x

molecules in subvolume n evolves according to

dP(x

,t)

dt = ˜ D ∑

m∈neighbours

{x

P(x

,t) − x

P(x

,t)}

+ g(x

− i)P(x

− i,t) + h(x

+ j)P(x

+ j,t)

− {g(x

) + h(x

)} P(x

,t)

(2.7)

The top row of eq. 2.7 describes the diffusion events, with rate constant ˜ D = D/l

, where l is the side-length of a subvolume.

Since eq. 2.7 is “just” a master equation, as eq. 2.3, but with an extended

state, it can be stochastically simulated using the same method as in sec-

tion 2.1.1. This method will, however, not be efficient, since the algorithm

re-samples event times in each iteration. In the RDME case, this means re-

calculations in a large number of subvolumes where nothing has changed

since the last iteration. To increase the efficiency of stochastic simulations

of the RDME, the next subvolume method (NSM) was developed [10]. This

method utilises the specific structures of the RDME in order to minimise the

number of times that event-times, i.e. the time until the next reaction or diffu-

sion event, have to be sampled. A general purpose tool for efficient stochastic

simulations of the RDME using the NSM, MesoRD, was developed as de-

scribed in paper VI. MesoRD was used for the stochastic simulations of the

Min-protein dynamics in paper IV.

2.2.2 The macroscopic limit of the mesoscopic description

As in the non-spatially dependent description above it is instructive to com- pare the mesoscopic, RDME, description to a macroscopic, deterministic de- scription. Restricting the analysis to one spatial dimension and taking the av- erage of eq. 2.7 yields [6]

d hx

i

dt = hg(x

)i − hh(x

)i + D

l

(hx

i − 2 hx

i + hx

i) (2.8) where i ± 1 are the two neighbouring subvolumes of subvolume i. The last part of eq. 2.8 is the Laplacian space discretisation of the second spatial derivative of Fick’s second law of diffusion [6]. The first two parts are the mean change due to reaction as in eq. 2.3. To go in the macroscopic limit of continuous space and time, where

∂ x(r, t)

∂ t = g(x(r,t)) − h(x(r,t)) + D ∂

x(r,t)

∂ r

(2.9)

would require small enough subvolume sizes, l, compared to the total system size, and at the same time large enough number of molecules in each subvol- ume to remove the stochastic effects of reactions.

Simulations of the mean-field description in eq. 2.8 have been used for the deterministic comparisons in paper IV, using an extended version of MesoRD.

2.2.3 The mesoscopic description in the microscopic limit

The results from the RDME description have been shown to diverge when the discretisation sizes become small [11, 12]. The reason for this divergence is the restrictions imposed when defining the RDME. In the RDME, each sub- volume is assumed to be well-mixed such that reaction probabilities only de- pend on the number of molecules in the subvolume where the reaction occurs.

In the case of bi-molecular reactions this leads to a lower limit on the subvol- ume sizes [10]. To complete a bi-molecular association event, A + B → C, in the mesoscopic description, the above restriction implies that the molecules first have to diffuse in order to collide and then react upon collision. Similarly, for the dissociation event, C → A + B, the molecules first have to dissociate and then diffuse far enough to be well-mixed. This restricts the number of models that can be implemented using the canonical RDME description. A correction method for irreversible reactions has been presented earlier [12].

Although a step in the right direction, this method is not based on a finer level of detail, but rather aims at retaining the copy number distribution described by the homogeneous master equation (2.1.1).

In the remaining part of this section is described how to define microscop-

ically correct reaction probabilities when the subvolumes can no longer be

considered well-stirred (for details, see paper V). The presentation here is re-

stricted to three spatial dimensions (see paper V for the corresponding two dimensional cases).

Consider the problematic bi-molecular reactions A + B

C

k_d

o

(2.10)

with the initial condition that the C-molecule is placed at the centre of a spher- ical volume. In the microscopic framework the C-molecule dissociates with rate constant γ into A and B. After dissociation the A molecule is placed in the centre and the B-molecule is placed on the reaction radius, ρ, away from the centre. The B-molecule can either diffuse away as a Brownian motion particle, with rate D, or react with A, with rate constant k, when it is on the reaction radius.

In the corresponding mesoscopic framework, space is discretised into con- centric shells of width h, starting from the reaction radius, ρ. The C-molecule dissociates with the rate constant q

, placing the B-molecule in the inner-most shell, spanning [ρ, ρ + h]. The B-molecule makes diffusion jumps between the shells. It can react with the A-molecule, with rate constant q

if it is in the inner-most shell.

In the conventional RDME, q

is set to the macroscopic reaction constant k

[13]. However, when h becomes small, q

should approach the microscopic association rate constant k to be consistent with the microscopic description.

Hence, here it is clear that the mesoscopic rate constants q

and q

depend on the size of the discretisation of the RDME. By requiring that the mean time to reach equilibrium for the scheme in eq. 2.10 is the same in both the micro- scopic and mesoscopic descriptions it is possible

to derive expressions for how q

and q

depend on the size of the discretisation, h (see the supplemen- tary material of paper V for details).

q

(h) = k/ [1 + α(1 + β )(1 − 0.58β )] (2.11) Here, β = ρ/(h + ρ) measures the spatial discretisation and α = k/4πρD measures how diffusion limited the reaction is. The dissociation rate constant is found as q

(h) = q

(h)/K, where K is the equilibrium constant for the reaction in eq. 2.10.

In figure 2.1A, the conventional RDME, using constant macroscopic, reac- tion rates k

and k

, is compared to the RDME with our new discretisation dependent reaction rates q

and q

in eq. 2.11. The new discretisation depen- dent RDME show discretisation independent results which are in excellent agreement with the microscopic description (compare the black and the red symbols, with the black solid line in figure 2.1).

1at least if you are Otto Berg

1

0.01 0.1

1 0.1

0.01 10

correct ka (n=2) qa (n=2) ka (n=10) qa (n=10) ka (n=100) qa (n=100)

Time

Pbound

A B

10^-3 10^-2 10^-1 10⁰ 10¹ 0

0.2 0.4 0.6 0.8 1

Time

l=2 , q(7l3) l=3.3 , q(7l3) l=5 , q(7l3) l=2 , q(l3) l=3.3 , q(l3) l=5 , q(l3) correct

Figure 2.1: Relaxation of the system in eq. 2.10 towards equilibrium. In (A), con- centric discretisation with n number of shells is used. In (B), Cartesian discretisation, with lattice size, l, is used. The conventional RDME implementation (black symbols) is compared to our new implementation using rates defined in eq. 2.11 (grey symbols) and the correct microscopic solution (black solid lines)

The next step is to extend the results, such that the corrections in eq. 2.11 can be used also for discretisation in Cartesian coordinates. The main ques- tion in this adaptation is how to choose appropriate values of h. In the example above, an A-molecule could react with B-molecules in the inner-most shell, i.e. at a distance ρ + h from the centre. Using cubic subvolumes, the same A-molecule can react not only with the B-molecules in the same subvolume, but also with B-molecules in all the neighbouring subvolumes. Imagine, for example, that the A-molecule is close to one of the edges of the subvolume.

To account for this uncertainty in position within the subvolumes, h is chosen such that the reaction volume is the same whether spherical or Cartesian co- ordinates with inclusion of neighbouring subvolumes is used. I.e. h is chosen from the equation (4π/3)(ρ +h)

= 7l

, where l is the subvolume side-length.

Figure 2.1B shows that the corrections in eq. 2.11 give excellent agreement with the microscopic description also for cubic subvolumes.

In steady-state situations, a change of binding kinetics will also change the steady-levels of the system. These effects are exemplified using the steady state system

A + B

C

k_d

o

/0

C B

/0 A

/0

(2.12)

The result of stochastic simulations of the system in eq. 2.12 for the conven- tional RDME (dashed lines) and the newly derived, corrected RDME (solid lines) is shown in figure 2.2. Using the conventional, uncorrected, version of the RDME the steady-state levels in the simulations depend on the discretisa- tion.

In conclusion, I have shown that using discretisation dependent association

and dissociation rates, the RDME produces microscopically correct kinetics

Number of C molecules

0 0.5 Time1 1.5 2

0 20 40 60 80 100 120 140 160 180

qA(7 ³) ka =2

=2.5 3.33

=4

=5

=10

=

Figure 2.2: Relaxation of the system in eq. 2.12 towards steady state. Cartesian dis- cretisation, with lattice size, l, is used. The conventional RDME implementation (dashed lines) is compared to our new implementation using rates defined in eq. 2.11 (solid lines). This figure was generated using an extended version of MesoRD, which allows for reactions in between subvolumes.

down to a spatial resolution corresponding to the reaction radii of the reacting

molecules.

3. Antibiotics as enzyme inhibitors and growth rate modulators

Antibiotic drugs are often enzyme inhibitors and their targets are often funda- mental processes of bacterial cells [14]. Numerous examples exist of antibiotic drugs which bind to the bacterial ribosome, thereby inhibiting protein synthe- sis [15]. Other common targets include the DNA replication machinery, the RNA polymerase and the bacterial cell envelope [16]. In order to understand the fundamental mechanisms of a particular drug’s inhibitory action, in vitro experiments of enzyme function in the presence of the drug can serve as a useful tool [17]. However, in order to understand how an in vivo phenotype, such as the growth-rate, is connected to a certain inhibitor-mechanism, the in vitro experimental data have to be put into the context of a growing cell.

In the first part of this section I will present the result of a model that was built to make the connection between in vitro observations of drug target bind- ing affinities, and in vivo observations of cell growth. More specifically the model gives an understanding of the in vivo experimental observation of drug target resistance masking, where drug target binding affinity differences do not give rise to growth rate differences in the living bacterial cell.

The model in the first part of the section does not make specific assump- tions regarding the mechanisms of drug binding to its target. However, specific knowledge of drug inhibitory mechanisms is a crucial component to make good, predictive models. In the second part of this section, the attention is therefore turned to what can be inferred regarding the binding mechanisms of enzyme inhibitors from in vitro steady-state kinetic experiments.

3.1 Drug target resistance masking

3.1.1 Antibiotic drug efficiency

The effectiveness of an antibiotic drug depends primarily on two things: (i) how well the drug can reach its target, and (ii) how well the drug binds when it is close to the target. Considering drug effectiveness as only the growth inhibitory effect of a single bacterial cell with an intracellular drug target, there are a number of things that affect drug effectiveness:

• The rate by which the drug enters the cell over the cell membrane, i.e. the

drug specific permeability of the cell membrane [18, 19].

• How well the drug binds its drug target when inside the cell, i.e. the binding affinity of the antibiotic drug to its intracellular drug target.

• The rate by which the cell actively pumps the drug into the growth medium, i.e. the drug efflux pump efficiency [20, 21].

• How long time the drug stays active when inside the cell, i.e. the drug degradation efficiency of intracellular enzymes.

3.1.2 Drug resistance

Resistance to antibiotics may arise in either of the above categories. Mutations in the porins have been shown to give drug resistance [22], suggesting reduced drug inflow. Drug target resistance mutations, which affect the binding affin- ity of drug to target, have been characterised for a large number of drugs [16].

Resistance mutations affecting the efficiency of efflux pumps have been char- acterised and shown to be of great importance for multidrug resistance (MDR) [20]. The MDR of bacteria with drug efflux pump mutations stems from the substrate promiscuity of efflux pumps [20]. Another resistance mechanism, for example in resistance towards β -lactams, is drug modification [16].

In a paper by Lovmar et al. [23], the drug target affinity difference of the an- tibiotic erythromycin binding to wildtype and resistance mutated E. coli ribo- somes affected in the ribosomal protein L22, was characterised. As expected, ribosomes with the drug target resistance mutation have a lower affinity for erythromycin then wildtype, and E. coli cells carrying the ribosomal mutation are less susceptible to erythromycin treatment. However, unexpectedly, the drug target resistance mutation does not lead to an in vivo growth difference, compared to wildtype, when the efflux pumps for erythromycin have been ge- netically compromised. Thus, the effect of the drug target resistance mutation in the ribosome has been masked and does not display a growth rate difference for cells with deficient efflux pump systems. Similar effects of drug target re- sistance masking in efflux pump deficient strains have also been observed by others [24, 25, 26].

To identify the mechanism of drug target resistance masking and to identify the criteria under which it occurs, we constructed a model of cells growing in a drug containing medium.

3.1.3 The model

Our model, described in the cartoon in figure 3.1, is an extension of the model

presented in Elf et al [27]. The model is given by the following set of differ-

ential equations (see figure 3.1 for parameter definitions)



 



 

 dr

dt = r

− r

k

− a

r

k

− µr

+ µr

da

dt = a

k

− a

k

+ r

− r

k

− a

r

k

− µa

(3.1)

We have analysed this model at steady state, where all the intracellular con- centrations and the growth rate are unchanged in time. At steady state, the relation between the normalised antibiotic inflow over the cell membrane and the normalised cell growth rate is described by

˜j

= (1 − ˜ µ )

˜k

+ ˜ µ

˜k

+ ˜ µ

˜k

µ ˜

+ ˜ µ

!

(3.2)

The following dimensionless variables have been used. ˜j

= a

k

/(r

µ

) is the normalised drug inflow, ˜ µ = µ /µ

is the normalised growth rate, ˜k

= k

/µ

is the normalised dissociation rate constant, and ˜k

= k

/(µ

r

) is the normalised association rate constant.

In paper I we have analysed this model and found the growth rate to be bistable under certain conditions. We also found that drug target masking is connected to growth rate bistability, such that an understanding of the bista- bility phenomenon is critical for the understanding of drug target resistance masking.

af

ka

kd

af

rf rb

kactive

af

ae

ae

af

kpassive

efflux pump

kpassive

ce^ll eⁿvelope

= antibiotic = antibiotic target kdeg

Figure 3.1: Cartoon illustrating drug flows over the cell envelope, drug inactivation

and target binding. Drugs enter the cell by passive transport (rate constant k

) and

exit by passive transport (rate constant k

), active pumping (rate constant k

),

drug degradation in the cytoplasm (rate constant k

) and dilution by growth (rate µ,

equal to µ

in drug free medium). Drugs bind the drug-targets with association rate

constant k

, and are released from drug-targets with dissociation rate constant k

. a

is the external, a

the free and a

the total drug concentration in the cell. The total,

free and drug-bound target concentrations are r

, r

and r

, respectively. In eq. 3.1,

k

= k

and k

= k

+ k

+ k

.

3.1.4 Steady state growth rate bistability

The steady state growth rate may become bistable as seen in figure 3.2. Growth bistability is characterised by the existence of a fast and a slow steady state growth rate, ˜ µ , for the same drug inflow, ˜ j

[28]. The two growth rate regimes (solid grey lines in figure 3.2) are separated by an unstable steady-state growth rate (dashed grey line in figure 3.2), which is never attained in the cell [28].

We further characterised the conditions under which growth rate bistability occurs. The conditions are low membrane drug permeability, inefficient drug efflux pumps, and high affinity, highly abundant drug targets.

We have in paper I given a careful analysis of when, where, and how growth rate bistability occurs. Here I will limit the description to promote an intuitive understanding of why growth rate bistability occurs under the conditions given above. This is done by following the arrows in figure 3.2, starting at low drug inflow in the fast growth regime. In the fast growth regime when the drug permeability is low and the drug has a highly abundant, high affinity, target, virtually all drug molecules that diffuse into the cell bind the drug target. Since the target bound drug molecules cannot diffuse out over the cell membrane, the main route of intracellular drug concentration reduction is via cell growth.

Thus, fast growth rate maintains the fast growth rate regime. Since the growth rate is monotonically decreasing with increasing drug inflow there will be a dramatic shift in the growth rate when the cell can no longer uphold the fast growth regime (vertical arrow down from fast to slow growth rate regime).

After the dramatic shift the cell will find a new slow growth rate, from which it cannot immediately escape if the drug inflow is lowered, since fast growth was a requirement for the fast growth regime (arrow to the left in slow growth regime). When the growth rate is fast enough the cell can revert to the fast growth regime (vertical upward arrow).

10^-3 10^-2 10^-1 10⁰ 10^-5

10⁰

Normalized growth rate,

Normalized drug inflow, j_in

Figure 3.2: Growth rate, ˜ µ as a function of drug inflow, ˜ j

. The relation between

growth rate and inflow is given by eq. 3.2 and are shown as grey lines (solid and

dashed). Stable steady state growth rate regimes (solid lines) are separated by unsta-

ble steady state growth rates. Arrows pointing downwards and to the right indicate

how the growth rate is changing during an increased drug inflow starting in the fast

growth regime. Arrows pointing upwards and to the left indicate how the growth rate

is changing during decreased drug inflow starting in the slow growth regime.

3.1.5 Drug target resistance masking

We suggest that drug target resistance masking is linked to steady-state growth-rate bistability. To highlight how efflux pump efficiency is linked to both bistability and drug target resistance masking, the steady-state growth-rates for three different binding strengths of drug to intracellular targets, in two different cases of efflux efficiency are shown in figure 3.3. In the efflux deficient strain (thin lines) the steady-state growth-rate is bistable, and there is no difference in growth-rate in the fast growth regime of the three different drug target binding affinities (for example, compare the influx

˜j

needed to reduce the growth rate to half of the drug untreated case, i.e ˜ µ = 0.5). This means that the efflux deficient strain shows drug target resistance masking in the fast growth regime. However, when the efflux pumps are proficient, there is no bistability and there are large differences in the steady-state growth-rate in the fast growth regime. The conditions for drug target resistance masking is thus the same as for steady-state growth-rate bistability. The underlying reason for drug target resistance masking is as follows. When the efflux pumps are deficient the main route of decreasing the intracellular drug concentration is via dilution by growth. This dilution effects both bound and free intracellular drug molecules equally, meaning that there will be virtually no difference between wildtype and drug target resistance mutants as long as the free concentration is negligible compared to the total. In the efflux pump proficient cases, on the other hand, the main route of decreasing intracellular drug concentration is via efflux pumps. The efflux pumps operate only on the free drug concentration, which will be different between wildtype drug target and drug target resistant mutants.

Normalized drug inflow, ^jin

Normalized growth rate, ₁₀^-2 ₁₀⁰ ₁₀² ₁₀^410-6

10^-5 10^-4 10^-3 10^-2 10^-1

10⁰ 0.5

0.004

Figure 3.3: Steady state growth rate, ˜ µ , as function of drug influx, ˜ j

, into the bac-

terial cells for drug-efflux proficient (thick lines) and deficient (thin lines) cells for

target wild type (black) and target resistance mutated (blue, 100 x wt dissociation

rate constant and green, 10 000 x wt dissociation rate constant) strains. For the target

wild type strain we use ˜k

= 10

, ˜k

= 0.01 and we use ˜k

= 10 in the drug-efflux

deficient and ˜k

= 10

in the drug-efflux proficient background. Solid and dashed

lines are stable and unstable steady-states, respectively. The upper and lower dotted

horizontal lines indicate normalised growth rates of 0.5 and 0.004, respectively.

3.1.6 Conclusions and outlook

In conclusion we have suggested a mechanism which explains the observed phenomenon of drug target resistance masking and we have also characterised the necessary conditions for masking to occur.

One necessary feature is that drug target resistance masking requires slow membrane diffusion. We immediately got feedback suggesting that the mem- brane permeability of passive membrane diffusion for erythromycin is about a 1000-fold higher than the parameters used in our examples in the discussion of paper I. The argument was based on measurements of the water/octanol parti- tioning coefficient for erythromycin [29] and the standard method of estimat- ing permeabilities from the partitioning coefficient [30, 31]. These numbers suggest that erythromycin equilibrates over the cell membrane with a charac- teristic time of 0.1 seconds. In a simplistic view, this means that in the efflux pump deficient mutants, the free internal drug concentration, a

, and the ex- ternal drug concentration, a

, are the same at all times, since this membrane diffusion is substantially faster then the cell growth rate. There are, however, a number of different observations suggesting that the drug concentration is not equilibrating very quickly.

• There is no growth difference in the efflux deficient strains. An equilibrium over the membrane would imply that growth rate differences should always be observed. This is however the reason why the model was constructed in the first place and therefore somewhat of a circular argument.

• The affinity of erythromycin to ribosomes is “too high” compared to the effect on growth in the efflux deficient cases. The affinity estimates of ery- thromycin is around 10 nM [23, 32] but also the efflux deficient bacteria grow happily with 1 µM erythromycin in the medium [23]. This argument is however weakened by the fact that we neither know the in vivo intracel- lular drug concentrations, nor the in vivo drug target affinities

.

• The only in vivo measurements of the rate of drug diffusion over the cell envelope that I could find, suggests that membrane diffusion of another drug, tetracyclin, is much slower [33] than the estimate given by the octanol/water partitioning coefficient [34]. The octanol/water partitioning does, however, only describe the inner membrane permeability and the comparison is therefore dependent on the outer membrane permeability.

• Goldman and Scaglione [32] suggest that erythromycin may use the porin pathway to traverse the outer membrane. This would invalidate all esti- mates of membrane diffusion rates based on octanol/water partitioning.

These points are all rather weak, making it hard to confirm that membrane diffusion is neither slow, as we are suggesting, nor fast, based on octanol/water

1Again, this is a reason to why we do need a model, since both these estimates are model dependent.

partitioning. This discrepancy does, however, strongly suggest the importance of a closer look at drug membrane permeabilities in a real in vivo situation

.

Together with Tanel Tenson at Tartu University we have tried to produce in vivo experimental data that might strengthen or overthrow our model. For example, the model suggests that the steady state growth rate is bistable and that the characteristic time to reach the new steady state growth rate after drug treatment should be very long when the drug concentration is close to the bifurcation point

. True steady state growth rate requires keeping cell cultures growing exponentially for long times, which turned out to be cumbersome to achieve. There might also be other issues complicating the situation. The response to erythromycin is heterogeneous, such that there is a wide range of cell lengths after erythromycin treatment. Our current macroscopic model can naturally not account for such heterogeneity and it is unclear if a stochastic model (see section 2.1.1) could reproduce this behaviour.

In the long term perspective, our model may serve as a tool for designing drugs with a reduced rate of drug target resistance evolution. This is espe- cially interesting in the light of the increasing effort that goes into finding ef- flux pump inhibitors [35], since our model predicts that drug target resistance masking requires deficient, or inhibited, efflux pumps.

3.2 In vitro characterisation of enzyme inhibitors

Characterisation of the mechanisms of action of enzyme inhibitors are readily done using in vitro systems. One of the standard ways of analysing enzyme in- hibitors is to measure the steady-state rate of product formation [36]. Inhibitor binding to enzyme can be described by the following classical scheme

E + S + I

k_I1



ES + I

q₁

o

k_I2



kc //

E + P + I

EI + S

q_I1

O

q₂

ESI

o

q_I2

O

(3.3)

By assuming that the substrate and the inhibitor equilibrates with the enzyme on the time-scale of product formation, the steady-state rate of product forma- tion, v, as a function of constant substrate, [S], and inhibitor concentration, [I], is given by

v

e

= k

([I])[S]

[S] + K

([I]) , or e

v = 1

k

([I]) + K

k

([I]) 1

[S] (3.4)

2As always, the only way to make progress is by disproving the model hypothesis.

3This is the point in figure 3.2 where the fast growth rate regime disappears, i.e. at the downward-pointing vertical arrow

where

k

([I]) = k

 1 + [I]

K



(3.5) K

k

([I]) = k

k

+ q

 1 + [I]

K



(3.6) Here e

is the total enzyme concentration and K

= q

/k

and K

= q

/k

are the dissociation constants for inhibitor binding to free and substrate bound enzyme, respectively. Some of the standard categories of enzyme inhibitors are attained in limiting cases of eq. 3.4. The limit when K

→ ∞ gives com- petitive inhibition, where the inhibitor can only bind the free enzyme. The limit when K

→ ∞ gives uncompetitive inhibition, where the inhibitor can only bind the enzyme substrate complex. K

= K

gives non-competitive in- hibition, where the inhibitor binds equally well to free and substrate bound inhibitor.

Although less restrictive, steady-state solutions of the scheme in eq. 3.3 have been known for a long time [37, 38], they have proved to be less use- ful for the analysis of enzyme inhibitors. The simplicity of the description in eq. 3.4 have made it ubiquitous in textbooks. However, one cannot conclude that substrate and inhibitor are equilibrated by studying the shape of eqs. 3.4- 3.6 alone. I will show below that even a perfect fit of experimental data to eq. 3.4 is not enough to show that the equilibrium approximation holds.

In paper II we revisit the steady-state solution to the scheme in eq. 3.3. We derive the steady-state solution of the scheme in eq. 3.3 using a new set of parameters. These allow for simple interpretations of special cases of the gen- eral steady-state solution, such as the substrate equilibrated solution described above and the rapidly equilibrating inhibitor description described below. In the paper we also describe what to expect in terms of precision and accuracy from regression analysis of experimental data for the general steady-state case and some special cases.

Here I will use an example of rapidly equilibrating inhibitors to show that

the interpretation of steady state enzyme kinetics data in the presence of in-

hibitors is not as straight forward as it may seem from textbooks.

3.2.1 Rapidly equilibrating inhibitors

By assuming that only the inhibitor is in equilibrium with the enzyme, the inhibitor binding is described by the scheme

q2

E+S+I

q1

ES+I ESI EI+S

KI2

KI1

kc

E+P+I

(3.7) Here brackets indicate quickly equilibrating parts in comparison to the re- mainder of the scheme. The steady state rate of product formation, v, as a function of the substrate concentration, [S], is in this case given by the same Michaelis-Menten equation as in eq. 3.4 ([37] and supplementary material of paper II). Here, however, the inhibitor dependent k

and K

/k

parameters are given by

k

([I]) = k

 1 + [I]

K



(3.8)

K

k

([I]) = R

 1 + [I]

K

 

1 + Qβ [I]

K





1 + β [I]

K

 (3.9)

The following parameters have been introduced:

• Q = q

/(q

+ k

), which describes how far away from equilibrium the sub- strate binding is. Q → 1 implies q

 k

and that substrate is equilibrated with enzyme on the time-scale of product formation. Q → 1 gives the equi- librated case above. Q  1 implies that ES far below equilibrium with E+S.

• R = k

(1 − Q), which is the k

/K

parameter for the enzyme in absence of inhibitor.

• β = k

/k

, which describes how much faster the substrate binds to inhibitor-bound then to free enzyme.

The functional form of K

/k

in eq. 3.9 has implications for the interpreta- tion of experimental data of enzyme inhibitors (see paper II). Three interesting cases appear:

• When β = 1, eq. 3.9 has the same functional form as eq. 3.6. This implies

that using eqs. 3.4-3.6 to analyse the data of a rapidly equilibrating in-

hibitors will give very high precision in the estimated parameters, but very

low accuracy since it may greatly miss-interpret the size of K

. In practise

this means that if substrate equilibrium is assumed, incorrectly, the size of

K

will be overestimated by a factor 1/Q. If Q  1 the estimated disso-

ciation constant, K

, becomes very large, and a non-competitive inhibitor will be incorrectly categorised as uncompetitive.

• When β = 1/Q the number of parameters in eq. 3.9 is reduced. Due to over- parametrisation there is low precision in parameter estimates, if eqs. 3.8 and 3.9 are used to analyse the experimental data.

• When β > 1 there is an ambiguity in how to assign the parameters in the model. The underlying reason for this ambiguity is as follows: when doing regression of eq. 3.9 on the form y([I]) = (1 + p

[I])(1 + p

[I])/(1 + p

[I]), it is not possible to separate the case where p

= β Q/K

and p

= 1/K

from the case where p

= 1/K

and p

= β Q/K

. In practise this means that a regression software may converge to any of two solutions depend- ing on the initial conditions. When β < 1, in contrast, β Q/K

is always smaller then 1/K

and the two cases can be separated.

In paper II we are suggesting a simple solution. All of the three problems

above may be resolved by measuring K

in a separate experiment. This mea-

surement can serve as a simple check to separate the cases of fast inhibitors

from equilibrated cases as described in the first point above. Furthermore,

beforehand knowledge of K

gives better parameter estimates in the second

point above. Finally, knowledge of K

can also resolve the parameter ambi-

guity described in last point above.

4. The effect of stochastic fluctuations on bacterial physiology and growth rate

Chemical reactions are stochastic, discrete events [3]. In the intracellular en- vironments of cells, such as E. coli, where cellular components exist in low copy numbers [39] the discrete nature of chemical reactions is intuitively ex- pected to generate large stochastic fluctuations in the amounts of intracellular components. This part of the thesis will investigate how these stochastic fluc- tuations of the intracellular components affect the physiology of the E. coli cell in three different cases.

In the first section (4.1) I will demonstrate how stochastic fluctuations in the flows of parallel anabolic biosynthetic pathways generate a fitness cost for the cell. I will show how this fitness cost can be partially remedied by the introduction of sensitive feedback loops and optimisation of the expression levels of the enzymes in the biosynthetic pathways.

In the two last sections I will describe how stochastic fluctuations and spa- tial distributions of proteins can affect the phenotype of cells. In section 4.2 the Min-system in E. coli is studied. The Min-system positions the cell division apparatus to the middle of the E. coli cell. In wildtype cells the Min-proteins oscillate from pole to pole. However, in a number of different mutants the os- cillation pattern is altered. We show that a stochastic reaction-diffusion model of the Min-system can be used to understand these mutant phenotypes.

In the third and last section (4.3) the transcription elongation of the RNA polymerases on ribosomal operons is studied. Using a DNA sequence depen- dent model for the RNA polymerase translocation rate two experimental ob- servations of the ribosomal operon in E. coli are explained.

4.1 Stochastic fluctuations in parallel anabolic pathways and their effects on cell fitness

To grow as fast as possible in any environment, an E. coli cell needs to opti-

mise its metabolism for maximal growth rate. For example, it has lately been

shown that uptake rates for various carbon sources together with oxygen in

E. coli conform to an in silico reconstructed metabolic network, where the

fluxes have been tuned for maximal biomass production per carbon input [40].

Furthermore, E. coli cells showing non-maximised biomass production on a specific carbon source have been shown to evolve towards solutions that max- imises biomass production [41]. Another view of E. coli growth optimisation is that, in order to grow at a maximal rate, the cell must distribute its resources according to a cost-benefit trade-off [42]. This type off trade-off has been ex- perimentally shown for the expression levels of the lac-operon in E. coli [43].

Here E. coli cells evolve towards an optimal cost-benefit trade-off for different concentrations of lactose.

For cells growing in amino acids lacking media, amino acid biosynthesis is fundamentally linked to cell growth rate. A production rate of amino acids below the maximal consumption rate by bulk ribosomes will define the rate of proteins synthesis, and thus the growth rate. The average translation rate has been shown to be governed by that amino acid production pathway, which has the smallest supply-to-demand ratio of amino acids [44], i.e. the limiting path- way. There are two reasons for why the translation rate depends on the flux of the limiting pathway alone. (i) The amino acids are produced in parallel path- ways and are sequentially incorporated into proteins, meaning that the fluxes of amino acids into proteins are stoichiometrically coupled to each other. (ii) The flux of amino acids into proteins saturate at high amino acid concentra- tions. This implies that amino acids in non-limiting pathways are incorporated into proteins at maximum rates, since they are supplied above their demand.

Thus, it is only the limiting pathway that affects the average translation rate.

Stochastic fluctuations in the production rates of amino acids show that the role of the limiting pathway is switching between the different pathways [44].

An open question is how stochastic fluctuations in the amino acids produc- tion pathways affect bulk protein synthesis and how average translation rate depends on the number of different amino acids that need to be synthesised.

The expression of amino acids biosynthetic enzymes is under feedback con- trol (e.g. [45, 46]). A second question is thus how the average translation rate is affected by these feedback control loops. Transcription initiation of the amino acid biosynthetic enzymes is controlled by repressors and ribosome depen- dent attenuation. Repressors regulate transcription initiation according to the intracellular concentration of the controlled amino acid [45]. The attenuation control system senses the translation rate of codons that are coding for the con- trolled amino acid [46]. Both of these control systems are negative feedback loops such that the transcription initiation rate increases when the supply-to- demand ratio of the amino acid decreases and vice versa. Negative feedback loops have been shown to decrease stochastic fluctuations in metabolite pools [47]. It has been suggested that one can understand feedback loops in terms of how well they manage to achieve optimal growth [48, 49, 50].

In the first part of this section of the thesis, I will describe how stochastic

fluctuations in parallel amino acids biosynthesis pathways affect the average

rate of protein synthesis. In the second part I will show how repressor con-

trol feedback loops can increase the average rate of protein synthesis and also

affect growth rate when stochastic fluctuations in parallel pathways are con- sidered. The details of this analysis is given in paper III and an overview is given below.

4.1.1 The model

We constructed a model of the amino acid turnover in a bacterial cell growing in amino acid lacking medium. The model is an extension of the model used in Elf & Ehrenberg [44]. The model contains the expression of amino acid biosynthetic enzymes, and the turnover of the enzymatically produced amino acids (fig. 4.1). Amino acid biosynthetic enzymes are expressed from their genes, via their mRNA. Each enzyme, in turn, produces one out of the twenty amino acids needed by an E. coli cell grown in amino acid lacking medium.

The amino acids are consumed on the ribosome, where they are incorporated into proteins. Before incorporation the amino acids need to be bound onto their corresponding tRNA. The transcription initiation of biosynthetic enzymes is regulated by repressor control systems sensing the intracellular amino acids pools. The model is implemented as a stochastic spatially homogeneous sys- tem (2.1.1).

E

AA1

mRNA1

E

RIB

AAN

Figure 4.1: Production and consumption of amino acids. When the bacterial cell is grown in amino acids lacking media, amino acids need to be produced by amino acid biosynthetic enzymes, E. All of the amino acids are produced in parallel to meet the demand set by the amino acid usage in proteins. Transcription of amino acid biosyn- thetic enzymes is regulated by negative feedback control systems.

4.1.2 Stochastic fluctuations and the rate limiting pathways

Before analysing the complete model shown in figure 4.1 (see section 4.1.3), I

will describe the underlying reasons for how stochastic fluctuations in parallel

amino acids biosynthetic pathways affect the average translation rate of the cell (see paper III and its supplementary material for details). To do this, I as- sume that all pathways are identical. This means that the catalytic rate of each amino-acid-producing enzyme is the same, the amino acids are incorporated with the same frequency into proteins, etc. In this view, the copy numbers of enzymes determine the production rate of amino acids in each pathway. This means that the enzyme pathway with the smallest copy number is limiting the translation rate, since it has the smallest supply-to-demand ratio. By further assuming that all enzymes are constitutively expressed, i.e. the transcriptional feedback control systems are temporarily removed from the model, and that expression of enzymes in each pathway is independent from the other path- ways, the stationary distribution of enzyme copy numbers in each pathway can be analytically derived [51]. By combining the known density function, f

, and distribution function, F

, of the enzyme copy numbers in each path- way

, the copy number of enzymes density function in the limiting pathway,

f

is given by

f

(E

) = N f

(E

)(1 − F

(E

))

(4.1) where N is the number of parallel pathways. The density function, f

, im- plies that one pathway has E

number of molecules, the other N −1 pathways have more then E

number of molecules and there are N different ways by which any one pathway can be rate limiting.

The average enzyme copy number in the limiting pathway, hE

i, decreases with increasing number of parallel pathways (figure 4.2A). This means that if the translation rate depends on the limiting pathway, the average translation rate will also decrease with increasing number of parallel pathways as seen in figure 4.2B. The cell growth rate µ = r hvi /ρ, is naturally proportional to both the average translation rate hvi, and the ribosome concentration, r, with a proportionality constant 1/ρ.

This means that for a constant given investment in ribosomes, a 1% decrease in the translation rate compared to the noise free case, gives a 1% decrease in growth rate, and thus a 1% increase in fitness cost due to noise. As seen in figure 4.2B the fitness cost due to noise increases with the number of noisy production pathways. The question which will be discussed in the following section is how the cell can reduce the noise induced fitness cost of parallel pathways.

4.1.3 Optimal growth rates

There are two ways by which the cell can reduce the fitness cost due to stochastic fluctuations in the amino acid production pathways. One way is to reduce the fluctuations in the enzyme copy numbers for each enzyme path-

1The enzyme copy number is distributed according to the Negative binomial distribution [51]

2ρ is the concentration of amino acids in protein [42]

<Elim>

Nr. of limiting pathways, N

1 10 100

600 700 800

900

A

Nr of limiting pathways, N^{0 10 20}

<v>/vmax

0.8 0.85 0.9 0.95

1

B

Figure 4.2: (A) The mean enzyme copy number, < E

> as a function of the number of pathways, N. (B) Normalised mean translation rate, < v > /v

, as a function of the number of pathways, N. v

is the macroscopically translation rate when each pathway is perfectly regulated and noise free. See paper III for details and parameter used.

way. To study the effect of noise reduction due to feedback control, I added repressor control systems back to the model. The repressor control was imple- mented with a simple Hill-type equation for the rate of transcription initiation, k

k = k

1 + ([AA] /K)

(4.2)

Here H is the Hill-coefficient, which determines the sensitivity of the repres- sor feedback, k

is the basal transcription rate in the absence of repressor, [AA]

is the amino acid concentration, and K defines the point of operation for the repressor control system. Eq. 4.2 implies that when the amino acid concentra- tion is high the transcription initiation rate is low, and vice versa.

Another way to reduce fitness-cost is to increase the investment in amino acid bio-synthesis pathways such that the ribosomes are saturated with amino acids also when there are stochastic fluctuations in the production pathways.

An increased investment in the enzyme pathway will, however, lead to a de-

crease in investment of all the other pathways, since the total protein invest-

ment in all pathways of the cell is constant. This restriction leads to a trade-off

between ribosome and enzyme system investments. I ran stochastic simula-

tions of the complete model in figure 4.1, which also included repressor con-

trol and limitations in the total resources available for the cell (see paper III

and its supplementary material for details). Figure 4.3 shows the cost-benefit

trade-off between enzyme systems and ribosomes, and how this is influenced

by repressor control. The constitutively expressed case (blue line in figure 4.3)

has its maximal growth rate at an over-expression of 25% of the enzyme sys-

tems compared to the investment which gives maximal growth rate in the noise

free, macroscopic, case (brown line in figure 4.3). This over-investment results

in a fitness-cost, where the maximal growth-rate is reduced due to reduced in-

vestment in ribosomes (compare the peak heights in figure 4.3). The repressor

control systems reduces this fitness-cost by reducing the noise in the enzyme