Improving landfill monitoring programs with the aid of geoelectrical - imaging techniques and geographical information systems Master’s Thesis in the Master Degree Programme, Civil Engineering KEVIN HINE Department of Civil and Environmental Engineering Division of GeoEngineering Engineering Geology Research Group CHALMERS UNIVERSITY OF TECHNOLOGY Göteborg, Sweden 2005 Master’s Thesis 2005:22 Mathematical modeling of bacterial populations under antibiotic selection pressure Master of Science Thesis in Engineering Mathematics and Computational Science Robert Andersson Department of Mathematical Sciences Division of Mathematics Chalmers University of Technology Gothenburg, Sweden 2013 MASTER OF SCIENCE THESIS IN ENGINEERING MATHEMATICS AND COMPUTATIONAL SCIENCE Mathematical modeling of bacterial populations under antibiotic selection pressure ROBERT ANDERSSON Department of Mathematical Sciences Division of Mathematics Chalmers University of Technology Gothenburg, Sweden 2013 Mathematical modeling of bacterial populations under antibiotic selection pressure ROBERT ANDERSSON c©ROBERT ANDERSSON, 2013 Department of Mathematical Sciences Division of Mathematics Chalmers University of Technology SE-412 96 Gothenburg Sweden Telephone: +46 (0)31-772 1000 Chalmers Reproservice Gothenburg, Sweden 2013 Mathematical modeling of bacterial populations under antibiotic selection pressure ROBERT ANDERSSON Department of Mathematical Sciences Chalmers University of Technology Abstract The development of bacterial resistance to antibiotic drugs have reached alarming proportions and is a serious threat to human health. As a consequence of improper antibiotic treatments, selection causes emergence of resistant bacteria. In order to avoid major medical issues, a better understanding of antibacterial dosage strate- gies is needed. In this thesis, the growth dynamics of bacterial populations under antibiotic treat- ment is modeled using coupled ordinary differential equations. Mutation dynamics is incorporated in the model as the bacteria have a certain probability to alter their genome in every new generation. The mutations affect the susceptibility of the bac- teria to the antibacterial drug. This leads to a selection pressure as new mutations develop. The aim was to find a good dosage strategy in terms of killing the bacteria with the least amount of drugs in the shortest time period possible. Here, three different antibacterial dosage strategies with constant concentration, a linearly decreasing concentration and exponentially decaying pulses were compared. Results suggest that the Minimum Inhibitory Concentration (MIC) of the slightly resistant bacteria acts as a threshold for what is a good dosage strategy. As long as the antibacterial concentration was kept over this level, no further resistant bacteria was allowed to develop and the system was destined to die. The strategy that most efficiently achieved this was the linearly decreasing regime. The model may be extended to usage within pharmacology in the future by incor- porating pharmacokinetics and pharmacodynamics. Keywords: bacteria, resistance, antibiotic, treatment, model, simulation i Acknowledgments I would like to express my greatest gratitude to my supervisors Assistant Professor Marija Cvijovic and Assistant Professor Erik Kristiansson for their enthusiastic involvement, patient guidance and assistance in my Master’s thesis. Without their help this thesis would not have been as successful nor fun as it has been. Gothenburg, August 15, 2013 ii Contents Abstract i Acknowledgments ii Contents iii Glossary and List of Abbreviations iv 1 Introduction 1 1.1 Antibacterial drugs and resistant bacteria . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Mathematical modeling in systems biology . . . . . . . . . . . . . . . . . . . . . . 3 1.3 Objective of thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 2 Model 5 2.1 Dynamical growth model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2.1.1 Mutations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.2 Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2.1 Ciprofloxacin - MIC and fitness . . . . . . . . . . . . . . . . . . . . . . . . 8 2.3 Simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 3 Results 10 3.1 Parameter space . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 3.2 Dosage strategy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 3.3 Simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 3.3.1 Constant concentration . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.3.2 Linear concentration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 3.3.3 Exponential decay . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3.4 Sensitivity analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 4 Discussion 23 4.1 Features of the dynamical model . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 4.2 Lack of precision in MIC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 4.3 Robustness of parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 4.4 Similarities between strategies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4.5 What is a good strategy? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 iii CONTENTS 4.6 Shortcomings of the model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 4.7 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 Bibliography 30 A Code 31 B Strategies 38 iv Glossary and List of Abbreviations Cmax is the highest antibacterial drug concentration [14]. ii, 5 F is the relative bacterial fitness in comparison to the wild type fitness which is 1. ii Nmax is the carrying capacity of the system. I.e. the maximal number of bacteria that the system can contain. ii kG is the growth coefficient of the bacteria. ii kmax is the maximum killing potential of the antibacterial drug. ii, 24, 25 p is the mutation rate of the bacteria per site and generation. ii, 24 AUC Area Under Concentration-time curve [9]. ii Bernoulli distribution is a one-parameter distribution Bern(p) that takes only two values, 0 and 1, with probability of 1− p and p respectively [10]. ii, 7 Binomially distribution is a two-parameter distribution Bin(n,p) which con- sists of n independent Bernoulli trials with probability parameter p[10]. ii, 7 CFU Colony-Forming Units. Measurement of bacterial number [17]. ii Ciprofloxacin is a type of fluoroquinolone [12]. ii DNA gyrase is a bacterial enzyme that prevents supercoiling in DNA replication [15]. ii, 2 v Glossary and List of Abbreviations E. coli (Escherichia coli) is a Gram-negative bacteria found in the intestinal tract of humans and other warm-blooded animals. Some strains of E. coli can cause sicknesses like urinary infections or food poisoning [15]. ii, 2 EC50 is the 50% effective concentration. I.e. the concentration in which 50% of the drugs maximal potential is attained [9]. ii Fluoroquinolones are antibacterial drugs that binds to and inhibits DNA gyrase and Toposiomerase IV under DNA synthesis [12]. ii Levofloxacin is a type of fluoroquinolone [17]. ii, 19 MIC (Minimum Inhibitory Concentration) is the lowest concentration that com- pletely inhibits visible growth of the organism as detected by the unaided eye after a 18- to 24-h incubation period with a standard inoculum of ap- proximately 105 cfu/ml [14]. ii, 1, 5, 10, 22–27 Normal distribution is a two-parameter distributionN (µ,σ2) and the most cen- tral in probability theory. According to the Central limit theorem it is the distribution of the sum of a large number of independent random variables. Examples of normally distributed variables can be the height of people, dis- tribution of IQ scores and so forth [10]. ii Pathogenic organisms are usually microorganisms that causes diseases [15]. ii, 2 PDI Pharmacokinetics/Pharmacodynamics Indices. ii, 5, Glossary: Pharmacoki- netics/Pharmacodynamics Pharmacokinetics/Pharmacodynamics Pharmacokinetics is in simplified terms what the body does to the drug while pharmacodynamics is what the drug does to the body [11]. ii Plasma is the liquid portion of the blood with cells removed and clotting proteins deactivated [15]. ii, 8 Topoisomerase IV is a bacterial enzyme that unlinks the two double-helical DNA strands from each other at DNA replication. It also prevents super- coiling, like DNA gyrase [15]. ii, 2 vi 1 Introduction T he development of bacterial resistance to antimicrobial agents is a serious threat to modern society and have recently reached alarming proportions [17]. Alterations in the genome of the bacteria is a cause of resistance, and the bacteria can e.g. acquire a genetic difference that affect its metabolic pathway. This may alter the cell wall and binding site of the antibacterial agent [6]. The problem of arising resistance can be correlated to suboptimal dosage strategies of antibiotics. This includes patients that not complete full treatments, usage of inappropriate dosages or wrong type of antibiotics [12]. All these issues lead to selection of resistant strains of bacteria. To overcome these problems and prevent unnecessary development of bacterial resistance, we need to improve our understanding of how to design good dosage regimens of antibiotics. An important tool is the usage of mathematical models. By modeling the growth dynamics of bacterial populations under antibacterial drug treatments, the effect of the drug can be analyzed at different time periods. This can give vital insight in design of dosage strategies. As a consequence, fewer experimental test and clinical trials are needed [17]. 1.1 Antibacterial drugs and resistant bacteria There are various different types of antibacterial agents. These can be categorized depending on their target area in the bacterial cell. The targets include inhibition of cell wall, protein and nucleic acid synthesis, interference with the bacterial membrane structure and disruption in the metabolic pathway of the cell [6]. On the other hand, random mutations in the bacteria can make it less susceptible 1 1.1. ANTIBACTERIAL DRUGS AND RESISTANT BACTERIA to the antibacterial agents. The mutated bacteria may e.g. acquire encoding enzymes that disrupts the antibacterial agents, alter the binding sites of the agents or generate efflux pumps that removes the agents from the cell. In this thesis, the main focus will be on a specific group of synthetic broad-spectrum antibiotics called fluoroquinolones and the drug ciprofloxacin in particular. The fluoroquinolones are used to treat a variety of different bacterial infections, like urinary tract infections, and they target the DNA synthesis of the bacteria [15]. When DNA is replicated, the bacterial enzyme DNA gyrase is used to remove supercoiling in the DNA along with the enzyme Topoisomerase IV, which also releases the the two replicated strands from each other. The fluoroquinolones dismantle these enzymes, which inhibits the replication. Escherichia coli (E. coli) is a gram-negative bacteria found in the intestinal tracts of humans and other warm-blooded animals [15] (see figure 1.1). Because of the bacterias substantial use in biochemical genetics, it was one of the main candidates to be the first organism with a fully sequenced genome [7]. Sequencing E. coli was a six years project which was finished in 1997. By then, the 4 649 221 base pairs of the DNA was sequenced and 4288 protein-coding genes were found. Pathogenic E. coli can cause infections in the urinary, pulmonary, enteric and nervous system, and was e.g. behind the outbreak of enterohemorrhagic E. coli (EHEC) in Germany between May and July 2011 [5]. Some infections caused by E. coli are treated with fluoroquinolones [7]. However, enzyme mutations located at the DNA gyrase and Topoisomerase IV in the E. coli decreases its susceptibility to the fluororquinolones. Figure 1.1: E. coli (Courtesy of Bioquell Inc). In this Master’s thesis, susceptibility data of ciprofloxacin will be used for wild type and mutated strains of E. coli. 2 1.2. MATHEMATICAL MODELING IN SYSTEMS BIOLOGY 1.2 Mathematical modeling in systems biology In systems biology, the strive is to analyze and understand cellular and organismal systems and structures [8]. A way to do this is to use mathematical models and methods that describe the given system or its properties. Mathematical models play a central role in systems biology. They can be categorized as deterministic and stochastic models. A deterministic system contains no form of randomness. Given a certain initial condition, the system will act the same and produce the same output every time. On the other hand, a stochastic system involve an amount of randomness and can give different output even though the input is the same. Biological systems are usually dynamical in their nature and are often modeled using ordinary differential equations. The principles of a modeling procedure using mathematics is illustrated in figure 1.2. Once the problem in question has been addressed, the modeling regime starts with a literature review to see what has been done in the field so far. With a schematic view of the system, key components can be highlighted and the first step of making a model of the system is done. Then, mathematical models are introduced to describe the dynamics of the system, often in the form of ordinary differential equations. Some of the parameters of the model may be unknown and for those estimates are needed. By using experimental data and e.g. goodness- of-fit, the parameters can be fitted and the model calibrated. With sensitivity analysis, robustness of the parameters can be determined to see in what magnitude different parameters influence the system. If a parameter is found to have little or no affect on the system, it may be set to a fixed number or less concern may be put on the accuracy of estimates of the given parameter [13]. To apply a mathematical model, an appropriate (depending on the specific model and taste) software is needed like e.g. matlab [19], Mathematica [21], Comsol [2] and so forth. Once the model can reproduce results from data, simulations with other setups can be done. In this way the model may be used to confirm, or reject, hypotheses or to predict features of the biological system. Of course, when working with models the pathway is not always this straight and one might need to go back and reevaluate some aspects of the model, use more input data for parameter fit and so forth. When it comes to mathematical modeling of antibacterial treatments, Bhagunde and colleagues [17] have developed a growth model for bacteria under antibacterial stress using ordinary differential equations. This model include an adaption func- tion describing the antibacterial resistance in a bacterial population. Even though the model is informative in the dynamics of the bacteria, it lacks dynamical infor- mation about individual subpopulations which consist of different mutated strains of the bacteria. It also lacks the stochastic nature of mutation development. 3 1.3. OBJECTIVE OF THESIS Figure 1.2: Schematic view of the constructing of a model [13]. 1.3 Objective of thesis The objective of this Master’s Thesis is to formulate a mathematical model that efficiently describes the growth dynamics of bacterial populations under antibiotic pressure. This will be done by applying modifications to the model by Bhagunde P. et al.[17] that describes effect of antimicrobials on heterogeneous bacterial pop- ulations. This model will in contrast to the model by Bhagunde P. et al. capture dynamics of several homogeneous subpopulations. Once the model is developed, a thoroughly exploration of its dynamics will be done in order to investigate the models properties and look at its sensitivity. With an antibacterial concentration function incorporated in the model, one can analyze different dosage strategies in order to find an efficient treatment to a bacterial infection. 4 2 Model T here are currently two main approaches for antibacterial pharmacody- namic/pharmacokinetic(PK/PD) modeling, those based on MIC analy- sis and those based on kill curves [14]. Depending on the given antibiotic agent used and its effect profile, different MIC models are used in order to determine a well-performing dosage strategy. Some antibacterial drugs show a time dependent killing and while they depend on the level of the dosage, they peak in effectiveness at a certain concentration and their effect is heavily depending on time they are present in the bacteria. Thus, the aim for such drugs is to have as long dosage time as possible. Other antibacterial drugs increase in effectiveness relative to the drug dosage and do not peak in effectiveness. In this case, the aim is to achieve the highest concentration possible instantly. The different MIC models used include Pharmacokinetics/Pharmacodynamics In- dices (PDI) [14][9] calculations like: • {t : C(t) > MIC}, the time where the concentration of the antibacterial drug is higher than the MIC. • Cmax/MIC, the fraction between the highest concentration and the MIC. • AUC/MIC, the fraction between the area under the dosage curve and the MIC. The other approach is to use time-kill curves. These models include the antibacte- rial concentration as a function of time and are aiming to capture the full dynamics of the bacterial populations. This approach is often more descriptive than the MIC- based one since it gives details about the effect of the antibiotics over time. In this project, a time-kill curve modeling approach is used together with MIC data 5 2.1. DYNAMICAL GROWTH MODEL for parameter fitting. The model describes the growth dynamics of eight bacterial populations. The main population consists of wild type E.coli while the other seven populations are mutations of the wild type. The mutations are located in the DNA gyrase (gyrA mutation), at two different sites, and in the Topoisomerase IV (parC mutation), at one site, giving a total of eight possible strains. The mu- tations affect the bacteria’s fitness, i.e. the bacteria’s potential to grow, as well as reduces its susceptibility to antibacterial drugs. In every new generation, a wild type E. coli can develop any of the three mutations, a bacteria with one mutation can develop any of the two remaining ones and so forth. A selective pressure will arise as a consequence of the fitnesses of the bacterial pop- ulations and their susceptibility to the antibacterial drugs. The more antibiotics that are added to the system, the greater the advantage will become to develop resistance. Given that the dosage of antibiotics will not wipe out the entire popu- lation, in the long run the wild type will be out-competed by mutations with lower susceptibility to the drug. 2.1 Dynamical growth model The growth dynamics of the bacterial populations are modeled by coupled ordinary differential equations. The i :th bacterial population, N i, follows a logistic growth given by dN i(t) dt = kGF iN i(t) ( 1− ∑ j N j(t) Nmax ) , for i = 1,...,8, (2.1) with kG being the growth rate of the bacteria, F i the fitness of the i :th population andNmax the carrying capacity of the system, i.e. the maximum number of bacteria that the system can hold. Without treating the populations with antibacterial drugs, the populations will be allowed to grow until they reach the maximum carrying capacity. Once reached, the populations with highest fitness will out- compete the other populations due to its higher growing ability. When applying an antibacterial drug treatment to the system, the bacteria start to die according to a three parameter Hill equation[18] given by dN i(t) dt = kmaxC(t) ECi 50 + C(t) N i(t), for i = 1,...,8 (2.2) with C(t) being the concentration of the antibacterial drug at time t. kmax is the maximal killing capacity of the drug and is attained when the drug concentration tend to infinity. ECi 50 is the antibacterial concentration in which 50% of the max- imal killing capacity is achieved. 6 2.1. DYNAMICAL GROWTH MODEL The kmax is considered to be drug specific and a measure of the drugs potency while ECi 50 is bacterial specific and a measure of the resistance of the i :th bacterial population. The full dynamics of the bacterial populations under a treatment of antibacterial drugs is then given by adding the growth of equation 2.1 and the the killing of equation 2.2 dN i(t) dt = kGF iN i(t) ( 1− ∑ j N j(t) Nmax ) − kmaxC(t) ECi 50 + C(t) N i(t), for i = 1,...,8. (2.3) Thus, all of the bacterial populations are modeled individually but their pro- gresses are dependent on each other due to the carrying capacity of the system. 2.1.1 Mutations In addition to the dynamics given by equation 2.3, the distribution of the bac- terial populations can be altered in every new generation due to development of mutations. If a mother cell lacks a certain mutation at a given site, the probabil- ity that its daughter attains that mutation follows a Bernoulli distribution with probability p. The mutation probability is considered to be fixed for all mutation sites. By summing over all daughters with mothers that lack a certain mutation, the amount which have attained that mutation in the next generation follows a Binomially distribution . With a generation time of approximately 20 minutes for E. coli [3], the mutation dynamics can be modeled by drawing a Binomial variable every 20th minute, corresponding how many daughters can gain a specific muta- tion. If e.g. there are N wild type E. coli at generation g, the number of daughter cells X that have developed a parC mutation in generation g + 1 is given by X ∼ Bin(N,p). (2.4) Similarly, the number of daughter cells X is computed for gyrA1 and gyrA2 mu- tations. Given that N is large enough, we can approximate X as a normal random variable with X ∼ N (Np,Np(1− p)). (2.5) The expected number of wild type bacteria that have developed e.g. a parC mu- tation in the n first generation, Y , is then given by Y = (N1 + ...+Nn)p (2.6) 7 2.2. DATA which is a sum over Normally distributed random variables and thus also Normally distributed. One can construct a 95% confidence interval for the expected number of bacteria that have developed a parC mutation in the n first generations Y ± zα/2√ n √√√√p(1− p) n∑ i=1 Ni, (2.7) where zα/2 is the α-quantile of the standard normal distribution. Here α = 0.05. 2.2 Data Traditionally in pharmacodynamics, the MIC has been used as a guideline when determining antibacterial drug dosages with the aim of having a Plasma drug con- centration that exceeds the MIC for as long as possible[4]. Recently, more detailed analyzes have been performed in order to get a better understanding of the dy- namics of the drugs. The MIC is, however, often used since it is easily measurable and as a consequence, there exist a lot of MIC data for different antibacterial drugs and bacteria. The MIC is by definition the antibacterial drug concentration in which the bac- terial population grows and dies at the same rate. In terms of equation 2.3, for population i and concentration C(t) = MIC this will imply: kGF iN i(t) ( 1− ∑ j N j(t) Nmax ) = kmaxMIC ECi 50 + MIC N i(t). (2.8) If we assume that the MIC is measured in conditions far from a saturated system, we would have that ∑ j N j(t) Nmax ≈ 0 and we can simplify the expression in equation 2.8 to kGF i = kmaxMIC ECi 50 + MIC . (2.9) From equation 2.9 we can estimate resistances EC50 for different bacteria given a certain drug and the MIC and fitness of the bacteria. 2.2.1 Ciprofloxacin - MIC and fitness In this study MIC and fitness values of the eight bacterial types for the antibacte- rial drug Ciprofloxcain, which is a type of fluoroquinolone, were used (table 2.1). 8 2.3. SIMULATIONS Bacterial strains include wild type (MG1655) and 7 different mutated strains. For example, strain LM378 have mutation S83L, meaning that the amino acid Serine (S) at the 83rd protein position is mutated into a Leucine (L). Table 2.1: Fitness and MIC of isogenic strains Strain gyrA1 gyrA2 parC MIC Fitness MG1655 - - - 0.016 1.00 LM378 S83L - - 0.38 1.00 LM534 - D87N - 0.25 0.99 LM792 - - S80I 0.016 0.99 LM625 S83L D87N - 0.38 0.97 LM862 S83L - S80I 1 0.98 LM1124 - D87N S80I 0.38 1.00 LM693 S83L D87N S80I 32 1.00 In the proposed model, the eight strains are divided into four populations depending on their susceptibility to the ciproflaxin. As can be seen in table 2.1, the LM792 strain shows the same susceptibility to the drug as MG1655, which is the wild type E. coli. These two strains make up the susceptible population. LM534 with a MIC of 0.25 and LM378, LM625 and LM1124, all with MIC of 0.38, will be considered as the slightly resistant population. The LM862 strain has a lower susceptibility and is the resistant population. Finally, LM693 has acquired mutations at all sites and makes up the heavily resistant population. 2.3 Simulations All simulations were carried out in matlab [19]. The system of coupled ordinary differential equations were solved using the ode45 solver. 9 3 Results T he following chapter contains the main results of this Master’s thesis. The result part will focus on the response in the dynamical model under different dosage strategies and levels of the antibacterial drug. Three dosage regimens will be presented and used under three different dosage levels. Finally, a sensitivity analysis will be performed on key parameters of the model. 3.1 Parameter space The data in table 2.1 was used to fit the parameters of the dynamical model. However, from equation 2.9, there are three unknowns for each data point (a data point consists of a MIC and fitness value). To reduce the parameter space to one unknown parameter per data point, assumptions were made on the parameter values of the growth rate kG and the maximal killing capacity of the drug kmax. The growth rate kG of the bacterial populations is considered to be constant over time and independent of the antibacterial drug concentration. Bhagunde et. al have fitted the growth rate for E. coli in a similar growth model to kG = 0.89 [17] and this is the value that will be used. The range of values for kmax was calibrated in the model such that the dynamics resemble real bacterial populations. The data was from the antibacterial drug ciprofloxacin and the kmax for this drug was set to 1.5. With these assumptions, one can fit bacterial resistances in terms of EC50 for all bacteria by using their MIC as following: 10 3.1. PARAMETER SPACE kGF i = kmaxMIC ECi 50 + MIC ⇐⇒ (3.1) ECi 50 = kmaxMIC− kGF iMIC kGF i , (3.2) for the i:th bacterial population. All parameter of the dynamical model given by equation 2.3 are presented in table 3.1 with the range of values that they attain. Parameter Description Values Units Assumption kG Growth rate of bacteria 0.89 cfu/h Same for all bacteria and all antibacterial dosages. kmax Maximal kill rate [1,2] h−1 Drug specific. EC50 50% of maximal kill rate concen- tration [0.0110, 21.9326] mg/L Bacterial specific. Measurement of the bacteria’s resistance of an antibacterial drug. Values fitted from MIC data. F Bacterial fitness [0.97, 1] Dimensionless Bacterial specific. The wild type has the highest fitness and a mutation can not increase fitness. MIC Minimum in- hibitory concen- tration [0.016,32] mg/L Nmax Carrying capac- ity of system 1010 cfu p Mutation rate [2·10−10,2·10−6] site / genera- tion Table 3.1: Parameter values of the dynamical model. 11 3.2. DOSAGE STRATEGY 3.2 Dosage strategy When setting up a treatment strategy for an antibacterial agent, the desired result is the killing of all bacteria and the prevention of any regrowth. In this sense a good strategy would be to simply treat with a very high dosage of antibiotics in a short time period killing the bacteria as fast as possible. However, severe side- effects can occur as a response to high dosage. Also, from an economical point of view, it is undesired to use an abundant amount of antibiotics. Thus, the key is to find a balance between these two criteria, i.e. killing all of bacteria should be achieved by lowest possible drug concentration at the shortest possible time. The dynamic of the system is modeled for different dosage regimens with the aim to find an optimal dosage strategy. The strategies are then compared in their ability to whip out the bacterial populations under the constraint that the AUC, i.e. the area under the concentration curve, is the same for the different drug dosage regimens. Optimally, one want to avoid any development of resistant bacteria. Once a colony of resistant bacteria has emerged, the antibacterial dosage needed to kill the populations is substantially increased. The response to the antibiotics, in terms of kill rate, for the different populations are given in figure 3.1 as a function of the antibacterial drug concentration. As a reference, the growth rate of the wild type is added. Once the kill rate is above this level, for a certain population, the population will die out over time. 0 0.5 1 1.5 2 0 0.5 1 1.5 Drug concentration (mg/L) K ill r at e Suceptible Slightly resistant Resistant Heavily resistant Growth rate of wild type Figure 3.1: Response to the antibiotics. Kill rate for the different bacterial popu- lations (y-axis) as a function of the antibacterial drug concentration (x-axis). The kill rates for the different populations will cross the growth rate line at the given population’s MIC. For the slightly resistant bacteria, MIC values are 0.25 and 0.38 mg/L respectively, for the subpopulations. This implies that as long as the antibacterial concentration is above 0.38 mg/L the slightly resistant bacteria will not grow and are kept at a controllable level. 12 3.2. DOSAGE STRATEGY The developed model was simulated for three different dosage strategies: (a) Con- stant antibacterial drug concentration (figure 3.2a), (b) linear decrease of con- centration from a high initial concentration (figure 3.2b) and (c) pulses with an exponential decay (figure 3.2c). 0 10 20 30 40 50 60 70 −0.5 0 0.5 1 1.5 Time (h) C on ce nt ra tio n (m g/ L) (a) Constant concentration 0 10 20 30 40 50 60 70 0.75 0.8 0.85 0.9 0.95 1 Time (h) C on ce nt ra tio n (m g/ L) (b) Linear decrease 0 10 20 30 40 50 60 70 0.4 0.5 0.6 0.7 0.8 0.9 1 C on ce nt ra tio n (m g/ L) Time (h) (c) Exponential decay Figure 3.2: Drug dosage strategies. Time(h) on x-axes and antibacterial drug concentration (mg/L) on y-axes for (a) constant (b) linear and (c) exponentially decreasing drug intake. The constant and linear strategies can be viewed as a drug treatment via intra- venous drip, where this could be the concentration profile in the blood stream. The exponential profile tries to mimic the absorption profile when taking antibiotics in the form of pills. 13 3.3. SIMULATIONS 3.3 Simulations The setup for the modeling of the dosage strategies consists of an initial population of 1010 wild type E. coli and no mutated strains. Mutations occur at a rate of 2 · 10−8 per site and generation. The simulations cover a time period of 300 hours. 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynamics 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population proportions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (c) Total population (d) Mutation dynamics Figure 3.3: Dynamics under constant antibiotic treatment. An example of how the simulations are illustrated. The subpopulations are represented with black, blue, green and red lines in (a), (b) and (d) for susceptible, slightly resistant, resistant and heavily resistant bacteria respectively. Time (h) on x-axes and on y-axis in (a) and (c) is the log10 of the population numbers, on y-axis in (b) is the proportions of the different subpopulations and in (d) is the number of mutated bacteria. To illustrate the dynamics of the system under antibacterial treatment, four different types of simulations are performed (example in figure 3.3). (a) popula- tion dynamics (figure 3.3a), illustrates the emergence of different populations over 14 3.3. SIMULATIONS time, (b) population proportions (figure 3.3b), depicts the proportions of the to- tal population that different population constitutes, (c) total population, (figure 3.3c) illustrates the total population number over time and (d) mutation dynam- ics, (figure 3.3d) shows the development of mutations over time along with a 95% confidence interval of the mutation numbers. 3.3.1 Constant concentration The constant concentration profile (figure 3.2a) is used at three different drug dosage levels: AUC = 75, corresponding to a constant drug concentration of 0.25 mg/l, AUC = 100, corresponding to 0.33 mg/l and AUC = 125, corresponding to 0.42 mg/l. The effect on the system of bacteria can be seen in figures 3.4 and 3.5. 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (b) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (c) Population dynam- ics (d) Mutation dynamics (e) Mutation dynamics (f) Mutation dynamics Figure 3.4: System response to treatment with constant antibacterial concentra- tion. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. As can be seen in figure 3.4, an AUC = 75 or AUC = 100 is not high enough to kill the bacteria while an AUC = 125 successfully kills all of the bacteria in approximately 175 hours. Heavily resistant bacteria are developed after approxi- mately 80 hours for AUC = 75 and after approximately 210 hours for AUC = 100. The curves of the susceptible, resistant and heavily resistant bacteria is very sim- 15 3.3. SIMULATIONS ilar when comparing figures 3.4a and 3.4b. The development of the resistant and heavily resistant bacteria is only occurring in a later stage of the treatment for AUC = 100. For AUC = 75 and AUC = 100, the slightly resistant and resis- tant bacteria are growing faster than they are being killed by the drug. However, they are out-competed when the heavily resistant bacteria are developed and are thereby eradicated from the system. In general, all dosages have quite similar ef- fect on all but the slightly resistant bacteria. When it comes to the mutation dynamics, the number of susceptible bacteria that develop slight resistance is approximately the same for all dosage levels. The slightly resistant and resistant bacteria develop more resistance for AUC = 75 than for AUC = 100 and for the AUC = 125 no slightly resistant bacteria develop resistance. 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (a) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (c) Population propor- tions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (d) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (e) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (f) Total population Figure 3.5: System response to treatment with constant antibacterial concentra- tion. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. In figure 3.5, the proportions of the populations as well as the total population are depicted. For a dosage level of AUC = 75 the system is saturated after approx- imately 100 hours (figure 3.5d) while for the AUC = 100 the system is saturated after approximately 200 hours. The AUC = 100 never led to a saturated system but is always declining until all bacteria are dead as can be seen in figure 3.5f. 16 3.3. SIMULATIONS From the proportion figures 3.5a, 3.5b and 3.5c can be seen that the subpopula- tions have different time periods when they dominate the system. For a dosage level of AUC = 75 and AUC = 100 the slightly resistant bacteria constitute the majority of the system until the resistant bacteria are developed and the system gets saturated (figures 3.5a, 3.5b, 3.5d and 3.5e). 3.3.2 Linear concentration As for the treatment with constant concentration, the linearly decreasing concen- tration (figure 3.2b) simulations are performed with AUC = 75, corresponding to an initial concentration of 0.5 mg/l that is dropping to 0 after 300 hours, AUC = 100, corresponding to an initial concentration of 0.66 mg/l and AUC = 125, corre- sponding to an initial concentration of 0.83. The dynamics can be seen in figures 3.6 and 3.7. 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (b) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (c) Population dynam- ics (d) Mutation dynamics (e) Mutation dynamics (f) Mutation dynamics Figure 3.6: System response to treatment with linearly decreasing antibacterial concentration. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. The linearly decreasing treatment kills all bacteria after approximately 40 hours for AUC = 100 and AUC = 125 as can be seen in figure 3.6b and 3.6c. An AUC = 75 is not enough to kill the bacteria (figure 3.6a). Figure 3.6a shows that resistant 17 3.3. SIMULATIONS bacteria is developed after approximately 180 hours and after approximately 210 hours the heavily resistant bacteria emerges. When it comes to mutation dynamics, the number of susceptible bacteria that develop slight resistant is almost the same for the different drug dosage levels (figures 3.6d, 3.6e and 3.6f). However, for an AUC = 75 there is also immense development of resistant and heavily resistant bacteria, figure 3.6d. 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (a) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (c) Population propor- tions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (d) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (e) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (f) Total population Figure 3.7: System response to treatment with linearly decreasing antibacterial concentration. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. Even though the AUC = 75 dosage level does not kill the entire bacterial population, the total population is very low (in the order of 10) within the time period of 50 and 100 hours, as can be seen in figure 3.7d. After approximately 200 hours the system with AUC = 75 is saturated. It is around that time period that system stops being dominated by the slightly resistant bacteria and bacteria with more resistance takes over (figure 3.7a). 18 3.3. SIMULATIONS 3.3.3 Exponential decay This treatment consists of injections every 8 hour (figure 3.2c). The concentration follows an exponentially decaying profile given by C(t) = Cmaxe −λτ (3.3) τ = tmod 8 (3.4) where Cmax is the maximum concentration that is attained just after an injection, λ is the decay constant and τ is the elapsed time since the last injection. The decay constant is related to the half-time t1/2 of the exponential decay as t1/2 = ln 2 λ (3.5) Nikolaou and colleges fitted an exponential decaying profile to data of Levofloxacin pharmacokinetics [17]. The best-fit value of the half-time was then 5.8 hours and this values is used in the simulations of exponentially decaying antibacterial con- centrations. Cmax is chosen to give one treatment with AUC = 75 (corresponding to Cmax = 0.39 mg/l), one with AUC = 100 (corresponding to Cmax = 0.52 mg/l) and one with AUC = 125 (corresponding to Cmax = 0.65 mg/l). The effect on the bacterial systems can be seen in figures 3.8 and B.6. 19 3.3. SIMULATIONS 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (b) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (c) Population dynam- ics (d) Mutation dynamics (e) Mutation dynamics (f) Mutation dynamics Figure 3.8: System response to treatment with exponentially decaying antibacterial concentration. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. The dosage levels of AUC = 75 and AUC = 100 are not killing the total bacterial population and heavily resistant bacteria are developed in both cases after approximately 80 and 160 hours respectively (figure 3.8a and 3.8b). For a dosage level of AUC = 125 the system is dying as can be seen in figure 3.8c. However, not on a 300 hours period, some slightly resistant bacteria are still alive at the end of the treatment. The population dynamics are similar when comparing the two lowest dosages with a shifting of the AUC = 100 in time in comparison to the AUC = 75. All but the susceptible bacteria grow faster than they die for these two dosage levels. However, with the arising of more resistant bacteria, a selection pressure is caused and eventually all but the heavily resistant bacteria are out-competed. 20 3.4. SENSITIVITY ANALYSIS As with the constant and linear dosage strategies, the mutation dynamics of the susceptible bacteria does not change significantly between different dosage levels. The slightly resistant and resistant bacteria develop more resistance for AUC = 75 than for AUC = 100 and does not develop for AUC = 125. 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (a) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (c) Population propor- tions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (d) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (e) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (f) Total population Figure 3.9: System response to treatment with exponentially decaying antibacterial concentration. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. As for the constant treatment (figure 3.5) the system gets saturated when the resistant bacteria takes over and dominates the system, depicted by figure 3.9a, 3.9d and 3.9b, 3.9e. 3.4 Sensitivity analysis When dealing with any mathematical model, sensitivity analysis is a vital tool in determining the robustness and structural weaknesses of the model. Sensitivity analysis can be divided into two classes, analysis concerning qualitative factors of the model and those concerning quantitative factors [20]. The focus here will be put on the quantitative factors which include e.g. perturbations in input variables and uncertainty in the model parameters. Uncertainty for the maximum killing 21 3.4. SENSITIVITY ANALYSIS capacity kmax and the mutation rate p was investigated using local methods, i.e. by keeping one of the parameter fixed while varying the other, in contrast to global methods in which all parameters are allowed to vary at the same time. 1 1.2 1.4 1.6 1.8 2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 k max (h−1) C on ce nt ra tio n (m g/ L) (a) System dependency of kmax −10 −9.5 −9 −8.5 −8 −7.5 −7 −6.5 −6 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 log 10 r C on ce nt ra tio n (m g/ L) (b) System dependency of p Figure 3.10: Sensitivity analysis. The threshold drug concentration (y-axes) needed to avoid heavily resistant bacteria in the long run for (a) different maxi- mal kill rates (x-axis) (b) different mutation rates (x-axis). In figure 3.10a the threshold constant concentration needed to avoid resistant bacteria in the long run, i.e. under a time period of 1000 hours for different values of kmax. Here the mutation rate is kept at a fixed level of p = 2 · 10−8. In figure 3.10b the kmax is kept at a fixed level of 1.5 while the mutation rate is varying. Determination of MIC often involves a 2-fold-dilution concentration series [4]. This means that starting with an initial sample with a certain amount of antibiotic drugs, the next sample will have a reduced concentration by a factor of two and the next after that a further reduction by an order of two and so on. The lowest concentration in the series that inhibits growth is then defined as the MIC. This method has an inherent level of uncertainty but its is practical and simple method of getting close to the real value [15]. If we let the C0 be the initial concentration and Ck be the k:th concentration then C0 = 2C1 = 4C2 = ... = 2kCk = .... (3.6) Thus, the potential error in the estimations of the MIC grows with the sample size by an order of 2. Assuming that the actual MIC is random number, the expected inherent error Ek for a drug with MIC = Ck is given by Ek = Ck − Ck−1 2 = C0(2 k−1 − 2k−2). (3.7) 22 4 Discussion I n this Master’s thesis, a dynamical model of bacterial growth by Bhagunde P. et al. has been modified in order to catch population behavior under selection pressure of antibiotic resistance. Parameters of the model has been fitted by MIC and fitness data of the antibacterial drug ciprofloxacin. In order to investigate the dynamics of the system and to find an optimal antibacterial dosage strategy, the effects of three different dosage regimens have been analyzed. Finally, the sensitivity of the model was analyzed in order to identify the robustness of the parameters. 4.1 Features of the dynamical model The model developed within this thesis incorporate both population dynamics of bacteria as well as genomics, in terms of mutation dynamics. This is a special model feature which there has been a lack of in the current field. The population dynamics are described by ordinary differential equations (equa- tion 2.3) where the bacteria follow a logistic growth and are killed according to a Hill equation once an antibacterial drug is added to the system. Mutations dynamics are incorporated as a stochastic factor that modifies the population dis- tributions at every new generation (described in section 2.1.1). The data used was obtained from a Etest (AB BIODISK, Solna, Sweden) [12] and the model relies heavily on the conditions at which this data was generated. Gen- erally, small changes in conditions may alter the results of biological experiments substantially. The saturation term in equation 2.3 enables selection in the system. Once the system reach a saturated state, the bacterial strain that most efficiently handles 23 4.2. LACK OF PRECISION IN MIC the antibacterial treatment will out-compete the other bacterial strains. This is shown in the proportion figures of all treatments where the antibiotic level is not enough to kill the system. Once the resistant bacteria start to develop, the slightly resistant bacteria will be out-competed and the system will be dominated by the resistant bacteria. This state is kept until heavily resistant bacteria develops. 4.2 Lack of precision in MIC The MIC is defined as ”the lowest concentration that completely inhibits visible growth of the organism as detected by the unaided eye after a 18- to 24-h incubation period with a standard inoculum of approximately 105 cfu/ml” [14]. The definition contains imprecise elements like detectability by ”the unaided eye” and a time period of ”18-h to 24-h”. As a consequence, different MICs from different labs or experiments can show essential variations. In the model used, the MIC plays an important role in parameter estimation (equation 3.1). As mentioned in section 3.4, when measuring the MIC one often use a 2-fold dilution series. The error level of this method can potentially be very big and thereby influencing the treatment strategy substantially. In this model, the MIC of the slightly resistant bacteria acts as threshold for what is a successful dosage level. Since the initial population of wild type E. coli is large in comparison to the mutation rate, slightly resistant bacteria will most likely develop at an early stage of the treatment. If the level of antibiotic drug is higher than the MIC of the slightly resistant bacteria, these bacteria will be destined to die and thus preventing further development into resistant bacteria. Thus, it is crucial to have an antibacterial drug concentration above the MIC of the slightly resistant bacteria to get a good treatment. Therefor the MIC of the slightly resistant bacteria need to be well-measured. 4.3 Robustness of parameters When simulating the population dynamics, assumptions where made on the two unspecified parameters kmax and the mutation rate p. In the sensitivity analysis in section 3.4, local methods was used to simulate the long-time behavior of the system in order to identify critical antibacterial dosage levels to avoid development of heavily resistant bacteria. When keeping the kmax fixed and varying the mu- tation rate p, the necessary antibacterial drug concentration needed exhibit small fluctuations for the different mutation rates (figure 3.10b). Thus, with respect to the long-time behavior, the system is very robust to fluctuations in the mutation rate. As long as the concentration is above the MIC of the slightly resistant bacte- 24 4.4. SIMILARITIES BETWEEN STRATEGIES ria, they will die and a very high mutation rate will be needed in order to develop resistance with significant probability. As a consequence of this, the model is not dependent on well-measured mutation rates to capture the dynamics of bacterial population. On the other hand, the kmax has naturally big influence of the system behavior (figure 3.10a). This is due to the central role the kmax plays in the killing of bacte- ria given by equation 2.2. The kmax decides what response a certain antibacterial drug concentration has on the system. 4.4 Similarities between strategies The exponential and constant strategy shows very similar behavior when it comes to the dynamics, as can be seen in figure 3.4 and 3.8. The constant strategy is more efficient as resistance is developed at a later stage of the treatment for AUC = 100 (figure 3.4b) and all bacteria are killed during the treatment for AUC = 125 (figure 3.4c) in contrary to the exponential treatments (figure 3.8b and 3.8c). This shows that the simulated drug has a time dependent killing profile, i.e. as long as we have a dosage level that is sufficient to kill the bacteria, an increase does not yield a substantially higher kill rate. However, as the exponential decaying profile fluctuates and gets under the MIC of the slightly resistant bacteria, this effect the killing. Even though the two dosage strategies show similar behavior, the exponential strategy suffers more in efficiency when it comes to its valleys in the drug profile than it gains from its peaks. Due to the similarities in dynamics between the two strategies, either one may be used depending on e.g. economical advantages. If a patient can take dosages by him or herself in the form of pills (corresponding to exponential dosage profile), this may be cheaper than having the patient hospitalized and on intravenous drip (corresponding to the constant dosage profile). 4.5 What is a good strategy? In the simulation section 3.3 the constant, linearly decreasing and exponentially decaying dosage profiles are tested for three dosage levels, AUC = 75, AUC = 100 and AUC = 125. The population dynamic results suggest that the linear decreasing dosage strategy is the most efficient one. As stated in the beginning of this thesis, the goal is to kill all of the bacteria with the lowest possible drug concentration at the shortest possible time. A key here is to kill the slightly resistant bacteria at an early stage of the treatment, before more resistant bacteria 25 4.6. SHORTCOMINGS OF THE MODEL are develops. To achieve this, a dosage level higher than the MIC of the slightly resistant bacteria is needed. And indeed, this is the case for all linearly decreasing profiles. However, with AUC = 75 the level drops under the MIC before the bacteria are killed. But with AUC = 100 and more, all bacteria are killed with the linear strategy while an AUC = 125 is needed for the constant and exponential strategies. According to the mathematical model, the best strategy is to have as high dosage as possible instantaneously to kill the slightly resistant bacteria before more resistance is developed. However, when dealing with real patients, other aspects need to be considered like the toxicity of the drug. A to high dosage may e.g. harm the liver of the patient and kill out other essential bacteria. 4.6 Shortcomings of the model ”As with all models, they are wrong but can still be helpful!” [16] A model can hardly ever capture the full dynamics of a real-life problem due to the immense complexity level of dynamical systems. However, a certain amount of features of the real life problem need to be included to make a model with predictive capa- bility. This is all a matter of balance between the complexity and the explanatory capability of the model. In this thesis, several simplifications of the real life problem has been made in order to start with a reasonably simple model: • The growth conditions are assumed to be constant over time meaning a system with ab libitum access to nutritions. A more realistic case may be to have a growth that is varying over time according to some stochastic process. • The mutation dynamics were modeled as discrete steps that occurred with a certain probability in every new generation. These mutations could then yield a lower susceptibility to the antibacterial drug. In reality, arise of mutations in the populations are more complex and there can be other factors leading to resistance beside mutations [6]. The parameter setting of the model is heavily dependent on reliable data. In order to understand the predictability of the model, experimental data is needed where growth dynamics could be compared to that of the model. A future aspect could be to incorporate pharmacodynamics and pharmacokinetics in the model. E.g, information of individual patients can in some way be included to construct a more patient based modeling profile. 26 4.7. CONCLUSIONS 4.7 Conclusions In this Master’s thesis, a mathematical model incorporating population dynam- ics and mutation genomics have been developed in order to investigate bacterial behavior under antibiotic treatment. This represents a first step in accurately de- scribing the dynamics of the populations, thus predicting development of heavily resistant bacteria, and suggesting a good antibacterial dosage strategy. The MIC describes the susceptibility of each strain to the antibacterial drug and is fundamental to the model. Slightly resistant bacteria will develop instantaneously and the key is to keep an antibacterial concentration higher than the MIC of these bacteria. In this sense, the linearly decreasing dosage strategy is the most efficient one used since it starts with a high dosage level and thus enables killing of the slightly resistant bacteria before less susceptible bacteria is developed. Because of the strong dependency of the measured MIC, the model is highly influenced by the experimental conditions in which the MIC was measured. To further evaluate the model, more experimental data is needed where population dynamics can be compared to that of the model. 27 Bibliography [1] Bernhard Ø. Palson. Systems Biology. Cambridge University Press, 2006. [2] COMSOL. Comsol multiphysics. [3] Edda Klipp. Modelling dynamic process in yeast. Wiley InterScience, 2007. [4] Elisabet I. Nielsen, Otto Cars, and Lena E. Friberg. Pharmacokinetic/phar- macodynamic (pk/pd) indices of antibiotics predicted by a semimechanistic pkpd model: a step toward model-based dose optimization. Antimicorbial Agents and Chemotherapy, 2004. [5] European Food Safety Authority (EFSA). E. coli: Rapid response in a crisis, 2011. [6] Fred C. Tenover. Mechanisms of antimicrobial resistance in bacteria. Ameri- can Journal of Infection Control, 2006. [7] Frederick R. Blattner, Guy Plunkett, Craig A. Bloch, Nicole T. Perna, Valerie Burland, Monica Riley, Julio Collado-Vides, Jeremy D. Glasner, Christopher K. Rode, George F. Mayhew, Jason Gregor, Nelson Wayne Davis, Heather A. Kirkpatrick, Michael A. Goeden, Debra J. Rose, Bob Mau, and Ying Shao. The complete genome sequence of escherichia coli k-12. Science, 1997. [8] Hiroaki Kitano. Systems biology: A brief overview. [9] Johan W. Mouton, Michael N. Dudley, Otto Cars, Hartmut Derendorf, and George L. Drusano. Standardization of pharmacokinetic/pharmacodynamic (pk/pd) terminology for anti-infective drugs: and update. Journal of Antimi- crobial Chemotherapy, 2005. [10] John A. Rice. Mathematical Statistics and Data Analysis. Cengage Learning, 2006. 29 BIBLIOGRAPHY [11] Leslie Z. Benet. Pharmacokinetics: Basic Principles and Its Use as a Tool in Drug Metabolism, pages 199–211. Raven Press, 1984. [12] Linda L. Marcusson, Niels Frimodt-Møller, and Diarmaid Hughes. Interplay in the selection of fluoroquinolone resistance and bacterial fitness. PLOS Pathogens, 2009. [13] Marija Cvijovic, Sergio Bordel, and Jens Nielsen. Mathematical models of cell factories: moving towards the core of industrial biotechnology. National Center for Biotechnology Information, 2010. [14] Markus Mueller, Amparo de la Peña, and Hartmut Derendorf. Issues in phar- macokinetics and pharmacodynamics of anti-infective agents: Kill curves ver- sus mic. Antimicrobial Agents and Chemotherapy, 2004. [15] Michael T. Madigan and John M. Martinko. Brock Biology of Microorganisms. Pearson Education, 2006. [16] Morris H. DeGroot. A Conversation with George Box. Statistical Science, 1987. [17] Pratik Bhagunde, Renu Singh, Kimberly R. Ledesma, Kai-Tai Chang, Michael Nikolaou, and Vincent H. Tam. Modelling biphasic killing of fluoroquinolones: guiding optimal dosing regimen design. Journal of Antimicrobial Chemother- apy, 2011. [18] Sylvain Goutelle, Michel Maurin, Florent Rougier, Xavier Barbaut, Laurent Bourguignon, Michel Ducher, and P Maire. The hill equation: a review of its capabilities in pharmacological modelling. Fundamental and Clinical Phar- macology, 2008. [19] M. The MathWorks Inc., Natick. Matlab 6.1. [20] Thomas Sumner. Sensitivity analysis in systems biology modelling and its application to a multi-scale model of blood glucose homeostasis. PhD thesis, University College London, 2010. [21] I. Wolfram Research. Mathematica. 30 A Code Here follows the MATLAB code that was used in this Master’s thesis: %main file close all %load data load data15 . mat data = data15 ; data ( : , 3 ) = 1 . 5 ; %time span T = 300 ; t = 0 ; %initial condition initial = zeros ( size ( data , 1 ) , 1 ) ; initial (1 ) = 10ˆ10; %matrices for output data tout = [ ] ; pout = [ ] ; %set the generationtime genTime = 1/3 ; %choose dosage strategy . dosage has three parameters for exponential strategy two for %linear and one for constant %exp − [ half time , dosage time , maximum dose ( c_max ) ] %lin − [ k , m ] %const − C strategy = ’ constant ’ ; dosage = 1/3 ; %flow of muations numberMut = zeros ( 8 , 1 ) ; expectedMut = zeros ( 8 , 1 ) ; varianceMut = zeros ( 8 , 1 ) ; %only call on mutation function if population size is greater than the %threshold mutationThres = 10ˆ4 ; 31 %index to store muation data index = 1 ; while t < T timeSpan = [ t t+genTime ] ; [ time population ] = ode45 ( @ (t , N ) odesys (t , N , strategy , dosage , data ) , timeSpan , initial ) ; timePoints = length ( time ) ; %temporary population and temporary statistics from mutations tempPop = zeros ( 8 , 2 ) ; tempTheta = zeros ( 8 , 2 ) ; %mutations for i = 1 : size ( population , 2 ) if population ( timePoints , i ) > mutationThres [ num theta ] = mutation ( population ( timePoints , i ) , i ) ; tempPop = tempPop + num ; tempTheta = tempTheta + theta ; end end %store mutation data if index ˜= 1 numberMut ( : , index ) = numberMut ( : , index−1) + tempPop ( : , 2 ) ; expectedMut ( : , index ) = expectedMut ( : , index−1) + tempTheta ( : , 1 ) ; varianceMut ( : , index ) = varianceMut ( : , index−1) + tempTheta ( : , 2 ) ; else numberMut ( : , index ) = tempPop ( : , 2 ) ; expectedMut ( : , index ) = tempTheta ( : , 1 ) ; varianceMut ( : , index ) = tempTheta ( : , 2 ) ; end %modify population according to mutations population ( timePoints , : ) = population ( timePoints , : ) + tempPop ( : , 1 ) ’ ; %store the data tout = [ tout ; time ( 2 : timePoints ) ] ; pout = [ pout ; population ( 2 : timePoints , : ) ] ; %set new initial conditions for next generation initial = population ( timePoints , : ) ; %if two low , population is considered dead for k = 1 : size ( population , 2 ) if initial ( k ) < 1 initial ( k ) = 0 ; end end %go to the next generation span t = t + genTime ; index = index + 1 ; end pout = abs ( pout ) ; %area under curve AUC ( strategy , dosage , T ) %population numbers in different time points pop96 = pout ( find ( tout >=96 ,1) , : ) ; pop168 = pout ( find ( tout >=168 ,1) , : ) ; pop300 = pout ( end , : ) ; %plot the data 32 plotData2 ( tout , pout , strategy , dosage , data , numberMut , expectedMut , varianceMut , T ) ; function dNdt = odesys (t , N , strategy , dosage , data ) %initialising dNdt = zeros ( length ( N ) , 1 ) ; fitness = data ( : , 2 ) ; K_max = data ( : , 3 ) ; EC50 = data ( : , 4 ) ; %coupled system of ODEs for i = 1 : size (N , 1 ) dNdt ( i ) = growthrate (t , N , fitness ( i ) , i ) − killrate (t , N , K_max ( i ) , . . . strategy , EC50 ( i ) , dosage , i ) ; end %growth rate of population , biofitness is considered to be constant between %wild types and mutations function growth = growthrate (t , N , fitness , target ) %growth rate constant , maximum population size Kg = 0 . 8 9 ; Nmax = 1∗10ˆ10; %this is wrong , it should be sum ( N ) growth = Kg ∗ fitness ∗ (1 − sum ( N )/ Nmax )∗ N ( target ) ; %kill rate of bacterial population function kill = killrate (t , N , K_k , strategy , EC50 , dosage , target ) %dosage strategy switch strategy case ’ constant ’ C = dosage ; case ’ pulse ’ C = pulseConcentration (t , dosage ) ; case ’ linear ’ C = linConcentration (t , dosage ) ; case ’ exponential ’ C = expConcentration (t , dosage ) ; end %sigmoidicity constant H = 1 ; %maximal kill rate , concentration to achieve 50% of maximal kill rate kill = CˆH ∗ K_k /( CˆH+EC50ˆH )∗ N ( target ) ; function [ out theta ] = mutation (N , target ) %output − first column represent population balance . second column %represent the cumulative mutation of the subpopulation out = zeros ( 8 , 2 ) ; %theta is the statistics , theta (1 ) = mu and theta (2 ) = sigmaˆ2 of the %normally distributed increments of the mutations theta = zeros ( 8 , 2 ) ; %mutation rate , 10ˆ−8 but scaled to make code faster %mutation factor M = 1 ; %mutation rate rate = M∗2∗10ˆ(−4); N = N / 10ˆ4 ; 33 %mutations . (1 −> 2 ,3 ,4 ) (2 −> 5 ,6) (3 −> 5 ,7) (4 −> 6 ,7) %(5 −> 8) (6 −> 8) (7 −> 8) if target == 1 out ( 2 , 1 ) = mutReturn (N , rate ) ; out ( 3 , 1 ) = mutReturn (N , rate ) ; out ( 4 , 1 ) = mutReturn (N , rate ) ; out ( 1 , 1 ) = out ( 1 , 1 ) − out ( 2 , 1 ) − out ( 3 , 1 ) − out ( 4 , 1 ) ; out ( 1 , 2 ) = out ( 2 , 1 ) + out ( 3 , 1 ) + out ( 4 , 1 ) ; theta ( 1 , 1 ) = 3∗N∗rate ; theta ( 1 , 2 ) = 3∗N∗rate∗(1−rate ) ; elseif target == 2 out ( 5 , 1 ) = mutReturn (N , rate ) ; out ( 6 , 1 ) = mutReturn (N , rate ) ; out ( 2 , 1 ) = out ( 2 , 1 ) − out ( 5 , 1 ) − out ( 6 , 1 ) ; out ( 2 , 2 ) = out ( 5 , 1 ) + out ( 6 , 1 ) ; theta ( 2 , 1 ) = 2∗N∗rate ; theta ( 2 , 2 ) = 2∗N∗rate∗(1−rate ) ; elseif target == 3 out ( 5 , 1 ) = mutReturn (N , rate ) ; out ( 7 , 1 ) = mutReturn (N , rate ) ; out ( 3 , 1 ) = out ( 3 , 1 ) − out ( 5 , 1 ) − out ( 7 , 1 ) ; out ( 3 , 2 ) = out ( 5 , 1 ) + out ( 7 , 1 ) ; theta ( 3 , 1 ) = 2∗N∗rate ; theta ( 3 , 2 ) = 2∗N∗rate∗(1−rate ) ; elseif target == 4 out ( 6 , 1 ) = mutReturn (N , rate ) ; out ( 7 , 1 ) = mutReturn (N , rate ) ; out ( 4 , 1 ) = out ( 4 , 1 ) − out ( 6 , 1 ) − out ( 7 , 1 ) ; out ( 4 , 2 ) = out ( 6 , 1 ) + out ( 7 , 1 ) ; theta ( 4 , 1 ) = 2∗N∗rate ; theta ( 4 , 2 ) = 2∗N∗rate∗(1−rate ) ; elseif target == 5 out ( 8 , 1 ) = mutReturn (N , rate ) ; out ( 5 , 1 ) = out ( 5 , 1 ) − out ( 8 , 1 ) ; out ( 5 , 2 ) = out ( 8 , 1 ) ; theta ( 5 , 1 ) = N∗rate ; theta ( 5 , 2 ) = N∗rate∗(1−rate ) ; elseif target == 6 out ( 8 , 1 ) = mutReturn (N , rate ) ; out ( 6 , 1 ) = out ( 6 , 1 ) − out ( 8 , 1 ) ; out ( 6 , 2 ) = out ( 8 , 1 ) ; theta ( 6 , 1 ) = N∗rate ; theta ( 6 , 2 ) = N∗rate∗(1−rate ) ; elseif target == 7 out ( 8 , 1 ) = mutReturn (N , rate ) ; out ( 7 , 1 ) = out ( 7 , 1 ) − out ( 8 , 1 ) ; out ( 7 , 2 ) = out ( 8 , 1 ) ; theta ( 7 , 1 ) = N∗rate ; theta ( 7 , 2 ) = N∗rate∗(1−rate ) ; end function out = mutReturn (N , rate ) out = ceil ( binornd ( round ( N ) , rate ) ) ; if isnan ( out ) out = 0 ; end %concentration shows a linear behavior function con = linConcentration (t , dosage ) %first input is the corresponding constant concentration 34 %second input is −1 for decreasing and 1 for increasing concentration k = 2∗ dosage (1)∗ dosage ( 2 ) /300 ; if dosage (2 ) == 1 m = 0 ; else m = dosage ( 1 )∗2 ; end %concentration con = k∗t + m ; if con <= 0 con = 0 ; end %concentration shows an exponential decrease function con = expConcentration (t , dosage ) %half time , maximum concentration at injection , time between doses tHalf = dosage ( 1 ) ; doseTime = dosage ( 2 ) ; cMax = dosage ( 3 ) ; %frequency , time since last dose lambda = log (2)/ tHalf ; elapsedDose = mod (t , doseTime ) ; %concentration con = cMax∗exp(−lambda∗elapsedDose ) ; function out = plotData2 ( tout , pout , strategy , dosage , data , numberMut , expectedMut , varianceMut , T ) figure ( ) wt = pout ( : , 1 ) + pout ( : , 4 ) ; res1 = pout ( : , 2 ) + pout ( : , 3 ) + . . . pout ( : , 5 ) + pout ( : , 7 ) ; res2 = pout ( : , 6 ) ; res3 = pout ( : , 8 ) ; plot ( tout , log10 ( wt ) , ’ k− ’ ,tout , log10 ( res1 ) . . . , ’ b− ’ ,tout , log10 ( res2 ) , ’ g− ’ ,tout , log10 ( res3 ) , ’ r− ’) legend ( ’ Susceptible ’ , ’ Slightly resistant ’ , ’ Resistant ’ , ’ Heavily resistant ’ , . . . ’ Location ’ , ’ Best ’ ) %legend ( ’ Susceptible ’ , ’ Slightly resistant ’ , ’ Location ’ , ’ Best ’ ) ylim ( [ 0 , 1 1 ] ) xlim ( [ 0 , max ( tout ) ] ) ; ylabel ( ’ log_ {10} N ( t ) ( cfu/ml ) ’ ) ; xlabel ( ’ Time ( h ) ’ ) %total population figure ( ) plot ( tout , log10 ( wt+res1+res2+res3 ) ) ylabel ( ’ log_ {10} N ( t ) ( cfu/ml ) ’ ) ; xlabel ( ’ Time ( h ) ’ ) ylim ( [ 0 , 1 1 ] ) xlim ( [ 0 , max ( tout ) ] ) ; %information of simulation %mutSum = sum ( mut , 1 ) ; %removed data of mutations , data instead showed in plot %mutText = [ ’ Mutations : ’ , num2str ( mutSum (1)+ mutSum ( 4 ) ) , ’ ’ , num2str ( mutSum ( 2 : 7 ) ) ) , ’ ’ , num2str ( mut ( 8 ) ) ] ; %k_max = [ ’ K_{max} = ’ , num2str ( data ( 1 , 3 ) ) ] ; %AUC = [ ’ AUC = ’ , num2str ( AUC ) ] ; %str = {k_max , AUC } ; 35 %ann = annotation ( ’ textbox ’ , [ 0 . 6 3 0 .6 0 .1 0 . 1 ] , ’ String ’ , str ) ; %set ( ann , ’ BackgroundColor ’ , [ 1 1 1 ] ) %fraction plot figure ( ) tot = wt + res1 + res2 + res3 ; plot ( tout , wt . / tot , ’ k− ’ ,tout , res1 . / tot . . . , ’ b− ’ ,tout , res2 . / tot , ’ g− ’ , tout , res3 . / tot , ’ r− ’) %legend ( ’ Susceptible ’ , ’ Slightly resistant ’ , ’ Heavily resistant ’ , . . . % ’ Location ’ , ’ NorthOutside ’ ) ylim ( [ 0 , 1 . 1 ] ) xlim ( [ 0 , max ( tout ) ] ) ; ylabel ( ’ Population proportions ’ ) ; xlabel ( ’ Time ( h ) ’ ) %mutation plot figure ( ) wtmut = numberMut ( 1 , : ) + numberMut ( 4 , : ) ; wtexp = expectedMut ( 1 , : ) + expectedMut ( 4 , : ) ; wtvar = varianceMut ( 1 , : ) + varianceMut ( 4 , : ) ; resmut1 = numberMut ( 2 , : ) + numberMut ( 3 , : ) + . . . numberMut ( 5 , : ) + numberMut ( 7 , : ) ; resexp1 = expectedMut ( 2 , : ) + expectedMut ( 3 , : ) + . . . expectedMut ( 5 , : ) + expectedMut ( 7 , : ) ; resvar1 = varianceMut ( 2 , : ) + varianceMut ( 3 , : ) + . . . varianceMut ( 5 , : ) + varianceMut ( 7 , : ) ; resmut2 = numberMut ( 6 , : ) ; resexp2 = expectedMut ( 6 , : ) ; resvar2 = varianceMut ( 6 , : ) ; %resmut3 = mut ( 8 , : ) ; gen = linspace (1 , round ( T ) , round ( T ∗ 3 ) ) ; hold on plot ( gen , wtmut , ’ k− ’ ,gen , resmut1 , . . . ’b− ’ ,gen , resmut2 , ’ g− ’ , . . . gen , wtexp , ’ k : ’ , gen , resexp1 , . . . ’ b : ’ , gen , resexp2 , ’ g : ’ ) ciplot (−2∗sqrt ( wtvar)+wtexp , 2∗ sqrt ( wtvar)+wtexp , gen , ’ k ’ ) alpha ( 0 . 1 ) ciplot (−2∗sqrt ( resvar2)+resexp2 , 2∗ sqrt ( resvar2)+resexp2 , gen , ’ g ’ ) alpha ( 0 . 1 ) ciplot (−2∗sqrt ( resvar1)+resexp1 , 2∗ sqrt ( resvar1)+resexp1 , gen ) alpha ( 0 . 1 ) hold off xlim ( [ 0 max ( gen ) ] ) ylim ( [ 0 8000 ] ) xlabel ( ’ Time ( h ) ’ ) ; ylabel ( ’ Number of mutations ’ ) %concentration plot figure ( ) switch strategy case ’ constant ’ C = @ ( t ) dosage ; fplot (C , [ 0 T ] ) case ’ exponential ’ C = @ ( t ) expConcentration (t , dosage ) ; fplot (C , [ 0 T ] ) case ’ linear ’ C = @ ( t ) linConcentration (t , dosage ) ; fplot (C , [ 0 T ] ) 36 case ’ pulse ’ C = @ ( t ) pulseConcentration (t , dosage ) ; fplot (C , [ 0 T ] ) end xlim ( [ 0 max ( tout ) ] ) xlabel ( ’ Time ( h ) ’ ) , ylabel ( ’ Concentration ( mg/L ) ’ ) ; function out = AUC ( strategy , dosage , T ) %calculate AUC switch strategy case ’ constant ’ out = dosage (1 ) ∗ T ; case ’ linear ’ out = dosage (1 ) ∗ T ; case ’ exponential ’ lambda = log (2)/ dosage ( 1 ) ; cmax = dosage ( 3 ) ; doseTime = dosage ( 2 ) ; out = ( T/doseTime )∗ cmax/lambda∗(1−exp(−lambda∗doseTime ) ) ; case ’ pulses ’ %something here end function ciplot ( lower , upper , x , colour ) ; % ciplot ( lower , upper ) % ciplot ( lower , upper , x ) % ciplot ( lower , upper , x , colour ) % % Plots a shaded region on a graph between specified lower and upper confidence intervals ( L and U ) . % l and u must be vectors of the same length . % Uses the ’ fill ’ function , not ’ area ’ . Therefore multiple shaded plots % can be overlayed without a problem . Make them transparent for total visibility . % x data can be specified , otherwise plots against index values . % colour can be specified ( eg ’k ’ ) . Defaults to blue . % Raymond Reynolds 24/11/06 if length ( lower)˜=length ( upper ) error ( ’ lower and upper vectors must be same length ’ ) end if nargin<4 colour=’b ’ ; end if nargin<3 x=1:length ( lower ) ; end % convert to row vectors so fliplr can work if find ( size ( x)==(max ( size ( x ))))<2 x=x ’ ; end if find ( size ( lower)==(max ( size ( lower ))))<2 lower=lower ’ ; end if find ( size ( upper)==(max ( size ( upper ))))<2 upper=upper ’ ; end fill ( [ x fliplr ( x ) ] , [ upper fliplr ( lower ) ] , colour ) S 37 B Strategies When investigating the effects of using the exponentially decaying dosage strategy in the Results chapter 3, antibiotics were injected every 8th hour. Here follows the results of using the same AUC with level 75, 100 and 125 but for injections every 6th, 12th and 24th hour. 38 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (b) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (c) Population dynam- ics (d) Mutation dynamics (e) Mutation dynamics (f) Mutation dynamics Figure B.1: System response to treatment with exponentially decaying antibacte- rial concentration. Injections every 6th hour. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. 39 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (a) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (c) Population propor- tions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (d) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (e) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (f) Total population Figure B.2: System response to treatment with exponentially decaying antibacte- rial concentration. Injections every 6th hour. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. 40 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (b) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (c) Population dynam- ics (d) Mutation dynamics (e) Mutation dynamics (f) Mutation dynamics Figure B.3: System response to treatment with exponentially decaying antibacte- rial concentration. Injections every 12th hour. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. 41 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (a) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (c) Population propor- tions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (d) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (e) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (f) Total population Figure B.4: System response to treatment with exponentially decaying antibacte- rial concentration. Injections every 12th hour. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. 42 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (a) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (b) Population dynam- ics 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 Time (h) lo g 10 N (t )( cf u/ m l) Susceptible Slightly resistant Resistant Heavily resistant (c) Population dynam- ics (d) Mutation dynamics (e) Mutation dynamics (f) Mutation dynamics Figure B.5: System response to treatment with exponentially decaying antibacte- rial concentration. Injections every 24th hour. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. 43 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (a) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (b) Population propor- tions 0 50 100 150 200 250 300 0 0.2 0.4 0.6 0.8 1 P op ul at io n pr op or tio ns Time (h) (c) Population propor- tions 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (d) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (e) Total population 0 50 100 150 200 250 300 0 1 2 3 4 5 6 7 8 9 10 11 lo g 10 N (t ) (c fu /m l) Time (h) (f) Total population Figure B.6: System response to treatment with exponentially decaying antibacte- rial concentration. Injections every 24th hour. In (a) and (d) the AUC = 75, in (b) and (e) the AUC = 100 and in (c) and (f) the AUC = 125. Susceptible, slightly resistant, resistant and heavily resistant bacteria are given by black, blue, green and red lines respectively. 44 Abstract Acknowledgments Contents Glossary and List of Abbreviations Introduction Antibacterial drugs and resistant bacteria Mathematical modeling in systems biology Objective of thesis Model Dynamical growth model Mutations Data Ciprofloxacin - MIC and fitness Simulations Results Parameter space Dosage strategy Simulations Constant concentration Linear concentration Exponential decay Sensitivity analysis Discussion Features of the dynamical model Lack of precision in MIC Robustness of parameters Similarities between strategies What is a good strategy? Shortcomings of the model Conclusions Bibliography Code Strategies