are.

a.

55%

ENE;

 

 

«ass...
.. 5.

..
5......
J15"...

aw.
.. 3
r fiti.‘
$¢.ahh.&n§
5%.?
3a:

a: 25....

2.3.“.
51‘ I.

 

v. 351. .
1.; 2. 55

 

THWS‘

7. (9’11

This is to certify that the

thesis entitled

Cellular Dynamics Simulation of Microbial
Chemotaxis in Porous Media

presented by

Jason R. Mondro

has been accepted towards fulfillment
of the requirements for

Master's degree in Chemical Engineering

 

 

V

 

I
Major professor

6/1/02

0-7639 MS U is an Affirmative Action/Equal Opportunity Institution

 

 

 

 

 

PLACE IN RETURN BOX to remove this checkout from your record.

LIBRARY

University

Michigan State

 

TO AVOID FINES return on or before date due.
MAY BE RECALLED with earlier due date if requested.

 

DATE DUE

DATE DUE

DATE DUE

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

6/01 c:/CIFIC/DateDue.p65-p. 15

 

fl,

CELLULAR DYNAMICS SIMULATION OF MICROBIAL CHEMOTAXIS IN
POROUS MEDIA

By

Jason R. Mondro

A THESIS
Submitted to
Michigan State University

in partial fulfillment of the requirements
for the degree of

MASTER OF SCIENCE
Department of Chemical Engineering and Materials Science

2002

ABSTRACT

CELLULAR DYNAMICS SIMULATION OF MICROBIAL CHEMOTAXIS IN
POROUS MEDIA
By

Jason R. Mondro

Chemotaxis refers to the directed migration of bacteria in response to a chemical
gradient. Engineered applications of chemotaxis offer the promise of solving difficulties
associated with in situ bioremediation for the removal of soil and groundwater pollutants.
Chemotaxis also plays an important role in the human immune system, and control of
chemotactic response may improve retention of surgically implanted medical devices.
Thus, an understanding of the mechanisms of chemotaxis, and strategies for its
manipulation is desired.

A cellular dynamics (CD) simulation method has been developed for the analysis
of chemotactic cell migration in porous media in response to gradients of consumable
chemoattractants. The simulation results are used to interpret experimental cell motility
measurements in heterogeneous porous aquifer solids. In CD simulations, macroscopic
cell population motility coefficients are calculated from the motion characteristics of
individual cells (i.e. swimming speed, tumbling frequency, turn angle distribution).

Using the CD simulation method, the motility and displacement of bacteria has
been investigated in bulk and porous environments, in response to consumable and
nonconsumable chemoattractants. The simulation results provide validation of several

experimental phenomena observed for migration of bacteria in porous aquifer solids.

DEDICATION

I would like to dedicate this thesis to my parents Thomas and Josephine for their support

and encouragement and especially to my wife Melissa for her patience and love.

ACKNOWLEDGMENTS

I would like to express my appreciation to Dr. Christian Lastoskie for all the help and
guidance in conducting this research. I would like to thank Dr. Mark Worden for his
insight into the experimental behavior of bacteria cells.

I would like to thank Dr. Syed Hashsharn for teaching me about the microbiology of
bacteria in a class that I really enjoyed and to Dr. Melissa Baumann for her inspiration
into the biomedical application of chemotaxis and implanted devices.

I would like to thank Ana Sarcheck for her help in running simulations and trying to keep
me organized. I would like to thank Sydney Forrester for his help in programming the
subroutine to place random particles in the model and for his expertise in using MatLab
to make the movies showing cell movement in front of the changing chemoattractant
concentration field.

I would also like to acknowledge Laura Booms and Caroline Roush for their

experimental photographs of interesting bacteria behavior.

TABLE OF CONTENTS

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

LIST OF TABLES ........................................................................... vii
LIST OF FIGURES ....................................................................... viii
NOMENCLATURE ......................................................................... xii
1. INTRODUCTION -- - - 1
1.1 BIOREMEDIATION BACKGROUND 1
1.1.1 Motivation and Severity of problem 1
1.1.2 Definition of problem--- 2
1.1.3 Relation to chemotaxis solution - - - 3

1.2 BIOMEDICAL BACKGROUND 5
1.2.1 Motivation and Severity of problem 5
1.2.2 Definition of problem 5
1.2.3 Relation to chemotaxis solution 6

1.3 OBJECTIVE -- 7
2. LITERATURE BACKGROUND .... -8
2.1 CHEMOTACTIC MECHANISM - 8
2.2 CURRENT EXPERIMENTAL PROTOCOLS -11
2.3 CURRENT MODELING PROTOCOLS - - - - 11
2.4 IMPROVED SIMULATIONS USING CELLULAR DYNAMICS 12
3. SIMULATION METHODOLOGY - - -- - - - - 15
3.1 CD SIMULATION MECHANICS - - - - - -15
3.2 GRID SPACING - - - -- 19
3.3 STABILITY 21
3.4 ASSUMPTIONS - - - - -24
3.5 DIMENSIONAerY 25
4. LINEAR GRADIENT NON-CONSUMABLE -- -- - 25
4.1 SINGLE CELL TRAJECTORY 26
4.2 EFFECTS OF VARYING CHEMOTACTIC SENSITIVITY COEFFICIENT - 28
4.3 EFFECTS OF VARYING CELL SWIMMING VELOCITY - - - -- - -32
4.4 EFFECTS OF VARYING THE TUMBLING PROBABILITY - - - 34
4.5 GAUSSIAN TURN ANGLE DISTRIBUTIONS - - - - - - -- 35
4.6 COMPARISON OF E. cou AND P. PUTIDA MOTILITY 37

 

4.7 EFFECTS OF THE MEAN AND VARIANCE OF THE TURN ANGLE DISTRIBUTION ..... 39

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

5. CONSUMABLE NON POROUS 43
5.1 INITIAL DISTRIBUTION 46
5.2 DIFFUSION 47
5.3 No DIFFUSION 48
5.4 ERROR ANALYSIS 50
5.5 ALTERNATIVE METHOD OF GRADIENT DETECTION 50
5.6 DEPENDENCE OF CONCENTRATION ON CHEMOTAxrs 55

6. MIGRATION IN POROUS MEDIA ............. .57
6.1 SLIT PORES--- 59
6.2 MULTISLIT BARRIER - - 64

6.2.1 Non-reflective barrier - - 68
6.2.2 Reflective barrier 76
6.3 HETEROGENEOUS 80
6.3.1 Migration into a strip of porous media 82
6.3.2 Migration in a bed of porous media 88

7. VISUALIZATION TOOLS ............. -- 97
7.1 EXCEL 97
7.2 ARRAY VISUALIZER - 98
7.3 MATLAB 98

8. CONCLUSIONS & FUTURE WORKS-- -- 99
8.1 CONCLUSIONS 99
8. 2 SUGGESTED TOPICS FOR FUTURE RESEARCH 102

8.2.1 3-D gradient detection and 3-D heterogeneous porous -- - 102
8.2.2 Multiple attractants/repellent consumption and or excretion ............ 102
8.2.3 Adsorption and energy tabulation to determine growth / death ........ 102
A PPE N DI X A ................................................................................ 1 03
APPENDIX B ................................................................................ 142
APPENDIX C ................................................................................ 143

REFERENCES .............................................................................. 145

Vi

LIST OF TABLES

Table 1. Values used in the simulations of a non-consumable linear gradient of a-
methylaspartate for E. coli. ............................................................................ 27

Table 2. Mean, median, and peak values of bacteria turn angle distributions. ................. 38

Table 3. Values used in the simulations starting with a uniform concentration of aspartate
consumed by E. coli. ...................................................................................... 44

Vii

LIST OF FIGURES

Figure 1. Schematic of flagella coordination by E. coli during Run and Tumble. ............ 9
Figure 2. Biasing of run lengths in the presence of a chemoattractant gradient. ............. 10
Figure 3. Simplified flow chart of model methodology. ................................................... 17
Figure 4. Schematic of the grid network with attractant being consumed based on cell
locations. ........................................................................................................ 18
Figure 5. Comparison of adjusting consumption based on grid spacing. ......................... 21
Figure 6. Cell trajectory comparison ................................................................................. 27
Figure 7. Snapshots of chemotactic and non-chemotactic cells in a linear gradient. ........ 28
Figure 8. Snapshots of cell location along a linear non-consumable gradient. ................. 30

Figure 9. Cell density profile in a linear gradient at varying chemotactic sensitivity

coefficients ..................................................................................................... 31
Figure 10. Chemotactic Sensitivity versus Motility and Mean Distance. ......................... 32
Figure 11. Velocity versus Motility and Mean Distance. ................................................. 33
Figure 12. Tumbling Probability versus Motility and Mean Distance .............................. 34
Figure 13. Normalized turn angle distribution varying mean. .......................................... 36
Figure 14. Normalized turn angle distribution varying variance. ..................................... 37
Figure 15. Turn angle distribution comparison of E. coli and P. putida ........................... 38
Figure 16. Cell density profile snapshot of E. coli and P. putida at 30 minutes ............... 39

Figure 17. Turn angle mean value versus motility and cell mean distance in a linear
nonconsumable gradient. ............................................................................... 40

Figure 18. Variance of the Turn Angle versus the motility and cell mean distance in a
linear non-consumable gradient. .................................................................... 41

Figure 19. Comparison of the mean turn angle in 2D, 3D against theoretical equation. .. 42

Figure 20. Comparison of the turn angle variance in 2D, 3D against theoretical equation.
........................................................................................................................ 42

Figure 21. Multiple Chemotactic waves of bacteria due to consumable chemoattractant.44

viii

Figure 22.
Figure 23.
Figure 24.
Figure 25.
Figure 26.
Figure 27.
Figure 28.
Figure 29.
Figure 30.
Figure 31.
Figure 32.
Figure 33.
Figure 34.
Figure 35.
Figure 36.
Figure 37.
Figure 38.

Figure 39.

Figure 40
Figure 41
Figure 42
Figure 43

Figure 44

Figure 45

Snapshots with Consumption Showing the cell density patterns. .................... 45
Radial cell density at different consumption rates. .......................................... 46
Effect of the consumption coefficient on motility with diffusion. .................. 48
Motility with diffusion and no diffusion on a logarithmic scale. .................... 49
Cell density profile where chemoattractant is depleted at center. ................... 51
Cell density profile where chemoattractant is not depleted at center. ............. 51
Comparison of different chemotactic sensitivities on the cell density profile. 52
Relative cell density profiles of gradient sensing methods with go, of 3.5E-5. 54
Relative cell density profiles of gradient sensing methods with x0 of 105E-4 55
Cell density profile and concentration profile due to consumption. ................ 56
The effect of concentration on cell density in a decaying step gradient. ......... 57
Schematic of cell trajectory through a bed of spheres. .................................... 58
Photograph of Pseudomonas KC migrating faster through sand core sample. 58
Simulation results of cells in a slit pore. .......................................................... 59
Effect of pore width on motility at different consumption rates ...................... 60
Motility at different population consumption rates. ........................................ 61
Effect of pore width on motility with different number of cells ...................... 62
Comparison of motilities in pores with consumption and no consumption. 63
Cell density profile for decaying step gradient in different pore widths. ........ 64
Migration into multiple 100-micron slit pores with a porosity of 33%. .......... 65
Migration into multiple lOO-micron slit pores with a porosity of 50%. .......... 65
Cell migration into multislits with a uniform initial concentration. ................ 69

Radial cell density profile and concentration profile with a uniform initial
concentration. ................................................................................................. 69

Cell migration into multislits with a lower initial concentration in the pores. 71

Figure 46. Radial cell density profile and concentration profile with a lower initial

concentration in the pores. ............................................................................. 71
Figure 47. Effect of chemoattractant concentration on penetration ratio. ......................... 72
Figure 48. Effect of chemoattractant concentration on the motility coefficient. .............. 72

Figure 49. Enhanced migration with initial lower concentration in multiple Slit pores... 73
Figure 50. Enhanced migration with initial uniform concentration in multiple Slit pores.74
Figure 51. Comparison of the penetration ratio on Short and long multislits. .................. 74
Figure 52. Comparison of the motility coefficient on short and long multislits. .............. 75
Figure 53. The effect of porosity on the penetration ratio and motility in multislits. ....... 75
Figure 54. Comparison of reflective and non-reflective walls of a ZOO-micron multislit. 76
Figure 55. Simulation results of reflective 900-micron multislit pores. ........................... 77

Figure 56. Simulation result of 100-micron pore width multislit with reflective walls.... 78

Figure 57. Comparison of 2D and 3D motion on penetration ratio. ................................. 79
Figure 58. Comparison of 2D and 3D motion on motility coefficient. ............................. 79
Figure 59. Cell trajectory through packed bed of 2D disks ............................................. 80
Figure 60. Cell trajectory through packed bed of cylinders ............................................. 81
Figure 61. Cell trajectory through a packed bed of spheres. ............................................. 81
Figure 62. Motility of different shape SOC-micron particles ............................................ 82
Figure 63. Simulation results showing cells migrating into porous region ....................... 83
Figure 64. Simulation result of cells starting out in a line along the porous region. ....... 84
Figure 65. Cell density profile for initial cell location of point versus line. .................... 85
Figure 66. Motility results for different particle diameter in strip. .................................. 85
Figure 67. Penetration ratio result for different particle diameter in strip. ....................... 86

Figure 68. Motility coefficient and Penetration ratio of 2D simulations with a porous
strip ................................................................................................................. 87

Figure 69. Migration in ordered 400-micron diameter 2D particles with 63% porosity. .88

Figure 70. Migration in random 400-micron diameter 2D particles with 64.6% porosity.
........................................................................................................................ 89

Figure 71. Migration in random 400-micron diameter 2D particles with 38.5% porosity.
........................................................................................................................ 90

Figure 72. Migration in random 400-micron diameter 3D particles with 51% porosity. 91
Figure 73. Effect of particle size on the motility coefficient in a bed of pores. ................ 92
Figure 74. Effect of porosity on the motility coefficient in 400-micron particle bed. ...... 93
Figure 75. Effect of porosity on motility coefficient in 500-micron particle bed. ............ 94
Figure 76. Cell density profiles for 2D and 3D motion with non-porous consumption. .. 96
Figure 77. Cell density profiles for 2D and 3D motion in pores with consumption ......... 97

Figure 78. MatLab output of tracking cell location and concentration in porous media. .99

xi

NOMENCLATURE

 

Number of cells

 

Time step

 

Length along positive X axis

 

Length along positive Y axis

 

Grid spacing

 

Tumbling probability (Berg and Brown 1972)

 

Cell swimming velocity (Strauss et a]. 1995)

 

Chemoattractant concentration

 

Dissociation constant (Rivero and Lauffenburger 1986)

 

m
U

N
o

Chemotactic sensitivity (Rivero and Lauffenburger 1986)

 

A
9-0
v

Mean displacement of the cells

 

N <"CJ ""
ragga .g<x 5

Random motility coefficient

 

Chemotactic motility coefficient

 

:1":

A
n

v

Position of cell at time t

 

Diffusion coefficient (Berg and Brown 1972)

 

Differential tumbling probability

 

366

Number of chemo receptor sites per cell

 

Turn angle resulting from a cell tumble

 

Direction vector of the cell

 

Consumption rate

 

Consumption rate

 

Length between nodal points

 

Weight distance of cell to nodal point

 

Variance

 

Turn angle

 

§¢Q€§W§MG

Mean turn angle

 

P
:2

Effective diffusion coefficient

 

m

Porosity or void fraction

 

 

c-l

 

Tortuosity

 

xii

 

1. INTRODUCTION

Two main areas for applications of chemotaxis include bioremediation of
contaminated water and soil and in the biomedical area such as in medical implants.
Chemotaxis has the potential to overcome several difficulties associated with in situ
bioremediation. By enhancing the dissemination of cells throughout a contaminated
aquifer, chemotaxis can help overcome slow advective transport and reduce fouling near
injection wells. By increasing microbial competitiveness, the chemotactic trait improves
the likelihood that an introduced microorganism will colonize the subsurface. Also, by
increasing cell motility, chemotaxis can extend the scope of the colonized region, thus
lowering capital costs and operating expenses by effectively reducing the number of
wells needed for delivery of cells and nutrients at the field site.

Chemotaxis also plays an important role in the recovery and successfulness of
surgery involving medical implants. Both bacteria and tissue cells utilize chemotaxis and
race to colonize surfaces of implants. The immune response and the healing process
produce chemotactic factors that are attractants for tissue cells and may act as repellants
to some bacteria. However, bacteria can migrate to the surface of an implant by different
methods. By modeling the response of motile bacteria in the body to chemoattractants
and chemorepellants, methods can be developed to repel bacteria from the implant site
thereby minimizing the chance of infection.

1.1 BIOREMEDIATION BACKGROUND
1.1.1 Motivation and Severity of problem

Over the years, solvents have other organic pollutants have been indiscriminately

dumped without consideration of the effects on the environment or on the people living in

the communities that are in close proximity to the dumpsites. According to a recent
report, “About 72 Million US. citizens live within 4 miles of a USEPA Super fund
National Priority List Site and 4.4 million people live within 1 mile.”(Penny and
Haldeman, 1994) We can no longer live by the mottos “out of sight out of mind” or
“dilution is the solution”. We have a moral obligation to dispose of chemicals properly
as well as to clean up our past mistakes. Therefore, how important is it to live in a clean
and safe environment? To what lengths should we go to achieve this goal?

It has been estimated that by using traditional decontamination methods, the cost
to “clean up” or remediate these sites could approach $1 trillion dollars (usgs. gov, 2000).
This estimate includes only the direct costs of clean up and not the costs associated with
the degradation of quality land, the associated medical problems, and the effects upon the
quality of life of local residents. Bioremediation, on the other hand, has the potential to
reduce costs, treat more sites, and remediate them faster than traditional methods. An
excellent example of a case that would benefit from bioremediation is contamination by
chlorinated solvents, because of the nature and toxicity of these solvents; priority in
remediation should be given to these chemicals.

1.1.2 Definition of problem

Bioremediation involves the use of biological methods to remove a pollutant from
an environment like soil or ground water. In this research we will be concerned with in
situ bioremediation, which involves the addition of non-native bacteria that have the
ability to break down the polluting contaminant. The inoculated bacteria need to compete
with indigenous bacteria for nutrients. A trait that can give it a competitive advantage is
chemotaxis. Chemotaxis helps overcome difficulty in bioremediation by enhancing

dispersion of bacteria into the contaminated area. This will help reduce the amount of

fouling and clogging at the injection point. It can also reduce the number of injection
wells and the overall cost.

An important aspect in bioremediation is the motility of the cells. If the cells are
non-motile or swim very slowly, the cells cannot reach the contaminated areas. A second
aspect is the cell diameter to the soil pore free space ratio. If the bacteria cell is too large
or the pores are clogged then the motility of the cell is hindered. Third assuming that the
cells can swim freely between the soil particles, how do we know which way they will
go? A bioremediation project will not be a success if the bacteria avoid the contaminant.
When choosing a bacteria cell for bioremediation a property of the cell that is needed, if
not required, is that the cell does not see the contaminant as a repellent.

1.1.3 Relation to chemotaxis solution

In bioremediation projects, a novel and effective pr0perty would be for the
bacteria to see the contaminant or a nutrient that is present in the contaminant as an
attractant. This property is called chemotaxis. With chemotaxis, bacteria will stay in the
areas of higher contaminant concentration for a longer period of time compared to a
bacteria utilizing only random motion (Duffy et al., 1995).

For some types of bacteria, oxygen is a chemoattractant, which can cause problems with
bioremediation (W idman et al., 1997). The oxygen concentration will be highest next to
the injection well and will drop off as the water containing the oxygen diffuses into the
soil. This gradient of oxygen will cause migration to the injection well, which can be in
the opposite direction of the increasing contaminant concentration. This mass migration
is also the cause of well plugging. This leads to the question “How do we know which
concentration gradient the bacteria will more likely pursue?” Studies of chemotaxis on

multiple concentration gradients Show the model that best fits experiments is a Simple

additivity of the gradients (Strauss, 1995). More recent research indicates that each
receptor doesn’t act independently and the activation of one receptor will increase the
activity of another receptor by lOO-fold and the whole array of receptors may be
involved. (Gestwicki, 2002).

Highly chlorinated chemicals in anoxic conditions undergo reductive
dechlorination, which involves removing chlorine atoms one at a time. Many bacteria in
anoxic conditions can degrade highly chlorinated chemicals, but cannot degrade lightly
chlorinated compounds. The problem with partial degradation is that the byproducts can
be more toxic than the original compound. An example is the degradation of TCE to
vinyl chloride. TCE (5 ppm) is a suspected carcinogen, while vinyl chloride (2 ppm) is a
known carcinogen (EPA, 1997). Another example of this is carbon tetrachloride, which
can be degraded to the suspected carcinogen chloroform (T atera et al., 1993). Even
though these lightly chlorinated chemicals can be degraded using aerobic bacteria, the
contaminated environment is not oxic so it won’t degrade. The addition of oxygen would
help to degrade the lightly chlorinated compounds, however the degradation of the highly
chlorinated compounds would stop.

Oxygen is a good electron acceptor and yields significant energy for cell growth,
metabolism and motility. One of the problems is that oxygen is not very soluble in water.
Only 9.1 mg/liter of oxygen is soluble in water at 20 °C (Newton, 1999). Since oxygen is
not very soluble it cannot be transported longer distances from the injection well. To
overcome the solubility problem more injection wells are needed with more additions of

oxygen added more often which will raise the cost of this remediation process.

Cells in an oxic environment grow faster than in an anoxic environments. The
logic, that bioremediation will be enhanced by a larger number of cells from faster
growth may have some faults. The faster rate of growth of cells has the possibility of
clogging the pore spaces between the soil particles, resulting in a sequestered cell
population that can’t get to the contaminant.

1.2 BIOMEDICAL BACKGROUND
1.2.1 Motivation and Severity of problem

Biomaterials and artificial prosthetics have been used by ancient Egyptians and
Greeks and were made out of plant material, wood and metals (Jacoby, 2001). Modern
devices including orthopedic joints, catheters, pacemakers, and cosmetic enhancements
are all made from a variety of different materials from metals to polymers and
transplanted tissues. A major cause of complication with implanted devices is infection,
which can result in sickness and death. While the infection rate in the 1960’s for people
undergoing a hip implant was 6.8%, it was still around 1% in the 1990’s. The infection
rate for knee and elbow implants is 2 to 9%, however it approaches 100% for indwelling
urethral catheters (Strickler and Mclean, 1995). The average cost for a knee prosthesis in
1991 dollars is $8,600 for a non-infected patient and $62,100 for an infected patient
(Strickler and Mclean, 1995). Since replacement of the implant is expensive and has
risks for the patient, more research needs to be focused on keeping the bacteria from
getting to the surface.

1.2.2 Definition of problem

Bacteria can get to the surface of implants through several different methods.
First, they can be initially present from improperly sterilized material or a breach in

sterile technique during operation. Second, airborne bacteria or indigenous bacteria from

the surrounding skin can enter the wound, even after the operation. Once bacteria invade
the body, they are carried by the flow of the tissue’s interstitial fluid and Brownian
motion. In this floating planktonic form the bacteria are more susceptible to antibiotics
and can be killed easier by the immune system (Gristina and Costerton, 1985). Bacteria
like Escherichia coli (E. coli) and Pseudomonas have flagella and actively swim to
implant surfaces. Motile strains of Streptococci have also been isolated contrary to the
belief that cocci bacteria are not motile (Van Der Drift et al., 1975). Once a thin layer of
bacteria or biofilm has formed on an implant surface, in most cases the only effective
treatment is removal of the implant (Chang and Merritt, 1994).

Bacteria alone in the body may not create an infection, but they become more
virulent once they adhere to a surface (Gristina et al., 1989). Bacteria in a biofilm can
grow faster and their metabolic activity can increase or decrease depending on the
nutrient availability (Hamilton, 1987). In this state, they are more resistant to antibiotics
and can have micro niches of multiple bacteria (Gristina et al., 1989). Over 99% of
bacterial biomass in nature exists as a biofilm, rather than free bacteria, indicating that
they would rather adhere to whatever solid surface they find (Gristina and Costerton,
1985). Once colonized, bacterial biofilms can spread and kill surrounding tissue. A
weakened immune system can increase the probability of infection so patients on
immuno-suppressant drugs are at high risk. Fortunately surfaces that are initially
colonized by healthy tissue tend to be resistant to bacterial infection (Ratner et al., 1996).

1.2.3 Relation to chemotaxis solution

Changes in extracellular pH from inflammation can affect the motility for bacteria
with flagella due to intracellular pH changes. High or low pH is a chemorepellant for

Escherichia coli and other bacteria. Around a pH of 3—4, flagella disintegrate into

subunits, while a pH between 6 and 8 doesn’t affect E. coli motility. Even at pH 7, when
in the presence of glucose Escherichia coli are not motile due to loss of flagella (Adler
and Templeton, 1967). The effect on motility of pH in the presence of different
attractants is not the same due to use of different receptors. The presence of different
chemicals can affect chemotaxis without binding to receptor proteins. Repellants,
membrane—permeant weak acids such as acetate and benzoate do not bind
chemoreceptors but shunt protons into the cytoplasm lowering intercellular pH affecting
chemotaxis (Slonczewski et al., 1982).

Heavy metal ions at very low concentration such as from metal implants inhibit
motility of E. coli, while amino acids act as chelating agents to stimulate motility (Adler
and Templeton, 1967). Unfortunately, build up of metal ions in surrounding tissue has an
undesirable cytotoxic effect and can lead to chronic inflammation and accelerated
corrosion of the implant leading to device failure. The presence of iron makes bacteria
resistant to killing by leukocytes indicating, its use as a biomaterial is not advisable
(Weinberg, 1974).

Both bacteria and bone cells utilize chemotaxis, thus methods should be exploited
to optimize the specificity in the design of implant devices to benefit tissue cells and
inhibit bacteria. For example, the porous Structure of hydroxyapatite (HA) can be used to
dispense antibiotics, chemoattractants for osteoblast, osteoclast and tissue cells, and
chemorepellants for bacteria.

1.3 OBJECTIVE

The objective of this research is to develop stochastic cellular dynamics (CD)
simulation calculations for cell motility coefficients from knowledge of the swimming

characteristics (swimming speed, tumbling probability, turn angle distribution) of

peritrichious microbes. Random motility simulations will be compared to theoretical
relations.

The effects of consumption on migration are determined and discussed in this
paper. Different mechanisms for gradient detection will be compared and analyzed. The
CD simulation will help determine the effects of porous media on the migration of
bacteria in a consumable chemoattractant. The boundary conditions involving the
interface with solid objects will be discussed and simulated to help explain observations
seen in experiments.

2. LITERATURE BACKGROUND

2.1 CHEMOTACTIC MECHANISM

Chemotaxis refers to the directed migration of organisms in response to chemical
concentration gradients. Chemotactic microorganisms migrate toward regions of high
chemoattractant concentrations at much faster rates than are realized through random
migration or convective transport alone (Witt et al., 1999). Chemotactic microorganisms
possess a significant ecological advantage over non-chemotactic competitors, in that
chemotactic Species are more readily able to secure nutrients needed for cell growth and
avoid potential toxins. Chemotaxis can therefore play an important role in engineered, in
situ, bioremediation systems for the treatment of soil/groundwater contaminant plumes.

Peritrichious bacteria cells such as E. coli and Salmonella typhimurium move in a
series of runs and tumbles. Bacteria swim by rotating their flagella clockwise, which
propels them forward. Bacteria change direction by reversing their rotation of their
flagella causing a tumble and resulting in a new swimming direction. Figure 1 shows the

flagella bundled during cell Swimming on the left and unbundled flagella during a tumble

on the right. This run and tumble motion can be simulated by a random walk (Berg and
Brown, 1972). The chemotactic effect is a result of the bacteria swimming for a longer
period of time up a perceived chemoattractant gradient before tumbling as shown in
Figure 2. The movement of bacteria unlike eukaryotic cells is actually a klinokinetic
response (Rivero et al., 1989; Ford, 1992). The overall directed motion of the bacteria
population results in a chemotactic effect. In chemoklinotactic simulations by Bombush,
the optimum turn angle that requires the least amount of energy to reach the attractant
source is between 40 and 100 degrees (Bombush and Conner, 1986). E. coli strain
AW405 has a mean turn angle of 68 degrees which is in the middle of the optimum range

(Berg and Brown, 1972).

 

Figure 1. Schematic of flagella coordination by E. coli during Run and Tumble.

Bacteria cells, which are on the order of 1 to 2 microns, are too small to sense a

chemical gradient along the length of its body so it can only bias its movement by

extending its run length. Bacteria use specialized proteins on their membranes called
receptors to measure the concentration of chemicals around them such as nutrients,
electron acceptors and toxins. Bacteria compare the concentration of a chemical in its
environment with the concentration sensed over the previous few seconds. As the cell
swims the temporal variations can be interpreted as a spatial gradient. The tumbling

probability is modulated according to the gradient detected.

 

Figure 2. Biasing of run lengths in the presence of a chemoattractant gradient.

Both aspartate and or-methylaspartate are chemoattractants for E. coli. Aspartate
can be metabolized and consumed by E. coli, while or-methylaspartate acts as an
attractant incapable of being metabolized. These two chemicals have similar Structures
enabling them to use the same receptor site on the cell. Even though their binding
disassociation constants (Kd) are different, they are the best candidates for chemotaxis
studies (T50 and Adler, 1974; Strauss et al., 1995; Mesibov et al., 1973). Furthermore,

aspartate is essential to biological systems and literature values of Kd have been reported.

2.2 CURRENT EXPERIMENTAL PROTOCOLS

Traditionally chemotaxis has been studied experimentally using capillary assays
and swarm plate experiments. Capillary assays are beneficial for screening bacteria to
determine if a particular chemical is an attractant or a repellant. The gradient at the
mouth of the capillary is difficult to characterize and control, thus making it difficult to
get accurate estimates of the cell motility and chemotactic sensitivity coefficients (Berg
and Turner, 1990; Ford, 1992; Marx, 2000; Vickers, 1981). Swarm plate experiments
provide a decent Visual representation of where and how quickly the bacteria migrate
from an initial position. In swarm plate experiments with sustained chemical gradients or
non-consumable (non-metabolizable) chemoattractants can’t easily be done, but this type
of experiment is amenable to characterize chemotaxis due to gradients produced by cell
metabolism. Visualization software is required for analysis of the cell density as a
function of time to extrapolate the motility parameters as well as fitting the mathematical
models. (Nikata et al., 1992; Strauss et al., 1995; Widman et al., 1997; Biondi et al.,
1998) An innovative method to analyze the behavior of bacteria in a fixed gradient has
been done using a Diffusion Gradient Chamber (DGC). In the DOC, constant
concentrations of chemoattractant is recirculated at one end and a buffer solution without
chemoattractant is recirculated at the opposite end of the chamber creating a fixed
gradient for studying consumable and non—consumable chemoattractants or
chemorepellants (W idman et al., 1997).

2.3 CURRENT MODELING PROTOCOLS

A number of mathematical models exist that have quantified motility parameters
and chemotactic sensitivity for Simplified physical conditions (Ford et al., 1991; Brosilow

et al., 1996; Rivero and Lauffenburger, 1986). One model uses a decaying step gradient,

11

which simulates a large l-Dimensional (1-D)chemoattractant gradient (Frymier et al.,
1994). Most of these models use different methods to solve different variations of the
cell balance equations, derived by Alt, to compute the unknown parameters. See Alt’s
publications for more details on the cell balance equations (1980). To further prove the
tractability of the CD method, Ford’s group used Monte Carlo method to solve Alt’s cell
balance equations (Frymier et al., 1993). Studies on consumable chemoattractants have
been done using cell balance equations or other methods, which don’t always quantify the
cell motion parameters (Marx and Aitken, 2000; Dillon et al., 1995; Widman et al.,
1997).

The motile behavior of a population of cells is characterized by both the mean
square distance the cell moves (motility) and the drift velocity (mean distance) along a
gradient. The motility is analogous to a diffusion coefficient indicating how dispersed
the cell population gets. The mean distance a cell population moves along a gradient is a
measure of its chemotactic response (Berg, 1983).

2.4 IMPROVED SIMULATIONS USING CELLULAR DYNAMICS

An alternative to solving cell balance equations is the cellular dynamics (CD)
method, which utilizes the cell motion parameters and independently calculates the new
positions of each cell. In addition, the algorithm considers the duration of time necessary
for a cell to tumble. The cell motion can be calculated using either two dimensions or
three dimensions, with the ability to calculate chemoattractant gradients in two
dimensions only. The CD model requires the chemotactic sensitivity as input. Therefore,
this parameter must be determined experimentally via a different method or taken from

reported literature values.

Three main independent scenarios were examined using the CD method. The first
case in focus is a fixed linear gradient with no consumption, which correlates with the
simplified chemotactic models previously done. This method is used to characterize how
the cell motion parameters affect chemotaxis. After first determining the optimal
chemotactic sensitivity to visually see displacement for the 30 minute simulation, several
of the cell motion parameters (i.e. cell velocity, mean turn angle, tumbling probability,
etc.) were adjusted to determine their effects on motility.

The second case modeled consumption starting with a uniform concentration of
chemoattractant. A significant advantage to the CD model is that the chemoattractant
concentration gradients can be independently measured in two dimensions. This is
required to model the chemotactic movement of bacteria that create gradients by
consumption.

The third area of investigation involved migration in porous media. Pore widths
and particle size were varied as well as the configuration of the porous region. The
model was adapted to simulate migration in Slit and heterogeneous pores with
consumption and diffusion of a chemoattractant and the simultaneous detection of the
gradient produced. Experiments with bacteria in both porous and non-porous media
Show patterns and other interesting phenomena. In an experiment with Pseudomonas KC
responding to nitrate, bacteria swam further through a porous sand core sample than
through the bulk phase (Figure 34). So we wanted to study migration in pores to
determine what was causing the accelerated migration. We also get some interesting

behavior without having pores. Other experiments Show migration away from the origin

in waves or hands of higher cell concentration resembling concentric rings. Both these
phenomena can be explained using our model for consumable chemoattractants.

The simulation model uses dimensionless variables in its calculations requiring
the input parameters to be scaled to a dimensionless number. The unit of length is sealed
with the value of 1 micron, which is the approximate length of an E. coli cell. The
scaling unit for time was selected to be 0.1 seconds; corresponding to the length of time a
bacterium tumbles (Berg and Turner, 1990). Concentration is scaled against the
dissociation constant to yield a dimensionless concentration. The consumption rate is
scaled with time from the input value of dimensionless concentration over time per cell.

Dillion uses a stochastic chemotactic model for the modeling of cells with
adsorption on to biofilms in porous media. These simulations track a few cells for only a
few seconds. Our CD simulations involve a large population of cells over 30 minutes and
longer giving longer-term dynamic effects of cell populations. Chemoattractant
concentration, individual cells and trajectories are discretely monitored in our model
rather than solving equations of motion. Dillion et a1. states that simplified chemotactic
models that ignore fluid effects give poor predictions of nutrient consumption (1995).
Nevertheless, our model is for non-flow systems and the chemoattractant does not have a
replenishing source or boundary to have significant fluid effects. We do however
account for diffusion of the chemoattractant. The hydrodynamic effect causing lower
concentration next to the cells can be accounted for by using a thin film mass transfer

coefficient with chemoattractant diffusing through to yield a surface concentration.

14

3. SIMULATION METHODOLOGY

3.1 CD SIMULATION MECHANICS

The bacteria runs are relatively straight and have a constant run time in the
absence of chemoattractant gradients (Brown and Berg, 1974). The basal tumbling
probability (p0) is the reciprocal of the straight line run time. When the cell is moving in
the direction of increasin g concentration the run lengths are extended. The basal
tumbling probability is adjusted based on the dot product of the direction vector (S) and
the concentration vector according to Equation 1 (Brown and Berg, 1974; Strauss et al.,

1995).

30
lo Kd

p 2p exp — SOVC (1)
0 v 2
(Kd +C)

 

E. coli only biases its tumbling probability when moving (up gradient) toward a
higher concentration of chemoattractant and returns to the basal tumbling probability
when moving (down gradient) away from a higher concentration of chemoattractant.
Salmonella typhimurium not only has a bias of longer runs when the concentration
increases, but exhibits shorter runs when swimming down gradient of the chemoattractant
as well (Macnab and Koshland, 1972).

The accepted expression for random motility with tumbles considered
instantaneous is defined by equation 2. The random motility coefficient 110 is
differentiated from the chemotactic motility coefficient by the subscript zero and is
calculated from a simulation of no gradient or consumption was compared tol .62E-6

cmz/s from the CD simulations. Using Equation 2 gives a random motility of 1.65E-6

15

cmzls resulting in a 1.85 % error compared to the random motility of 1.62E-6 cmzls
calculated by the CD simulation (Liu, 1996; Biondi et al., 1998). A smaller random
motility coefficient of 1.4E-6 cm2/s was calculated from a simulation using 0.05 second
time step and 0.1 second delay for a tumble. Neglecting the inclusion of the actual time
of a tumble results in a difference of 15.1 % in the random motility coefficient compared
to the actual random motility. Using a time step of 0.1 instead of 0.05 seconds yields a

decreased random motility with the addition of a 5% error.

_ v2 (2)
_ 2p0(l—<cos<¢>>>

 

#0

In equation 2, the coefficient 2 refers to the dimensionality of the system. The
average turn angle (i) that is calculated from the turn angle distribution is 792° in the
model. The motility coefficient is only invariant when the gradient is unidirectional and
non-changing or linear. Consumption from a single point yields a mean square
displacement slope that changes over time because the gradient also changes as the cells
thin out in a circular wave from their initial location. At large consumption rates where
cells get stranded behind, large deviations from a linear slope occur because part of the
cell population is moving with random motility and the other portion of cells respond to a
changing gradient as the cell wave front proceeds.

By using a cellular dynamics model, the stochastic nature allows the chemotactic
response of individual cells to be randomly biased giving a realistic population behavior
at larger population sizes. A simplified flow chart of the model methodology is shown in

Figure 3.

l6

A random number is generated to determine if the cell tumbles or runs. If the
random number is less than the new tumbling probability, the cell tumbles, thus staying at
its present location. A new direction is then picked based on its current direction and its
turn angle distribution. If the cell runs, the direction stays constant and a new location is
determined based on the cell velocity. This process is repeated for each cell for each time

step.

i

initial position
and direction

1‘

calculate new
tumbling probability

 

 

 

 

 

 

 

 

 

Yes No

 

 

choose new direction calculate new
based on turn angle position

probability dist. r1(t+At) = r1(t) + At

 

 

 

 

 

 

 

 

 

 

 

+5
I

Figure 3. Simplified flow chart of model methodology.

The motility coefficient (It) is obtained from the mean displacement of the cell

population (52(0). The cell displacement is averaged over the number of cells (N)

(Duffy, 1995).

1 N
Q“) =fi2

[de + dyz]
i=1

(3)

17

Einstein Relationship for 2 dimensions

Q(r) —+ 4p r
as r —> oo
(4)
Expression for the motility coefficient
1 ( am) (5)
l1 = — —
4 t

The concentration of chemoattractant is tracked at nodal points on a 2-

Dimensional grid with the origin located at the center (Figure 4).

 

Figure 4. Schematic of the grid network with attractant being consumed based on cell locations.
The four nodal points around the cell are updated for consumption of each cell
and are weighted (W) based upon its relative distance from each node. The consumption
is modeled as a constant or zeroth order reaction rate. A Michaelis-Menten kinetic rate
could be used as well for concentration dependent consumption of nutrients. The rate of
consumption in these simulations is assumed constant where the change in concentration
per time is defined as k.

—=k (6)

The chemoattractant concentration field is updated at the end of each time step for
diffusion. The edges of the grid can be chosen to be periodic or non-flux depending on
the details of the Simulation. Equation 7 shows how the concentration at each point is

updated due to cell consumption.

C(x, y,t + A!) : C(x, y,t) —W(x, y) k At

(7)
3.2 GRID SPACING

The 2-dimensional grid works well for keeping track of the chemoattractant
concentrations at fixed nodal points as Shown in Figure 4. Grid spacing is extremely
important in obtaining an accurate concentration profile. If the grid spacing is too small,
the cell will swim past several nodal points at each time step, creating gaps in
consumption along the cell path. This could lead to an erroneous concentration profile.
Alternatively, if the grid spacing is too large the nodal square for consumption will be
unrealistically large as well. As a result, the calculated gradient would be located far
from the cell position. Furthermore, cell movement will approach random motility at
high grid spacing values. An error analysis was completed on adjusting the grid spacing
ranging from values of 24 pm to 40 um. The error percentage, though small, is
noticeably higher for the simulations completed using a grid spacing of 40 um, and lower
for a grid spacing of 24 um.

The effects of the grid spacing on the chemotactic motility coefficient were
determined for Simulations with consumable chemoattractant. Grid spacing would
appear to have a Si gnificant effect on the gradients calculated using the nodal points.

When the total consumption on a concentration basis at the four closest nodal points is

19

kept constant and the grid spacing is changed, the shape of the motility trend is similar to
the case of changing consumption rates because the consumption rate is effectively
changed. The consumption rate needs to be constant per volume or area of the grid
space, so adjustment of the changing area of the nodal square needs to be incorporated to
conserve the rate of mass consumption. An area based on 20-micron grid spacing was
assumed for adjusting the consumption rate in 2 dimensions at different grid spacing
using the expression in equation 8. The units on the consumption rate k and k0 are
(dimensionless concentration/(number of cells * seconds)). The consumption rate R has a
constant mass basis, while kO has a constant concentration basis for use at the 4

surrounding nodal points.

M (8)

 

By adjusting the consumption rate based on a constant mass basis yields a
constant motility as the grid space changes as shown in Figure 5. Diffusion is
significantly larger than the consumption rate so the gradients smooth over the nodal
points making the simulations less dependent on grid space. This allows us to have more
flexibility in setting the size of the grid space when running simulations with

consumption.

20

 

 

 

 

 

 

 

 

2.5E—5

2.0E—5 ~ - L *
E
5 155-5 . -- —
8’
€1.0E-5~ - ,s- ~~~ AAAAA~~
e
2

5.0E—6 ~ -‘ n * -_- -— . +Mass bas'sk

-9— Concentration basis ko
0.0E+O T . .
O 10 20 3O 4O
Gridspace (pm)

Figure 5. Comparison or adjusting consumption based on grid spacing.
3.3 STABILITY

The criteria for numerical stability were determined using the center finite
differencing method to approximate diffusion. To increase the margin of stability, the
grid spacing should be as large as possible, while keeping the time step to a minimum.
With our algorithm, the time step was fixed to the length of time a bacteria takes to
tumble (approximately 0.1 seconds), leaving only the grid spacing available to adjust.
The grid spacing is used to determine the gradient that the bacterium senses. It is
believed that bacteria sense a gradient temporally rather than Spatially. Therefore, the
grid spacing must correlate with the distance the bacteria swims during the time it takes
to notice a change in the concentration of its surroundings (Macnab and Koshland, 1972;
Brown and Berg, 1974). As in the case of consumption, the grid spacing cannot be
arbitrarily set.

The minimum grid spacing required for numerical stability was calculated to be

21.1 um using a time step of 0.1 seconds. This value was calculated by solving a

21

concentration balance to determine the grid Size and time step that balance the diffusion
in and out of a nodal point. This nodal point has a lower concentration than the four
surrounding points. Using both the basal tumbling frequency and cell velocity, the
approximate average distance a cell travels before it senses a gradient is 18.3 um. Based
upon the calculations and an error analysis, a grid spacing of 24 um was chosen for our
Simulations to ensure stability as well as give a reasonable estimate of the cells ability to
sense gradients.

Our simulations are based on a 2-D grid with concentration saved at the nodal
points. The consumption at each nodal point is based on the number of cells present in
each grid space and how close they are to the nodal point. Diffusion is accounted for by
using an explicit finite difference scheme utilizing equations 9 and 10.

32c' 32c

C(r + At) =C(r) + DAt —— + — (9)
3x2 8Y2

 

 

82C Cj+1,k “ZCjJe +Cj—l,k
2 = 2 (10)
EX DC
A diffusion balance was done to determine the stability criteria for the grid
spacing and time step giving equation 11.
< _
Ccenter + dC center '— Csurrounding dC surrounding (11)

The change in concentration due to diffusion at the center point and the 4

surrounding points are defined by the following expressions.

22

dC = 4dC (12)

 

 

center surrounding
center D 2("surrounding — 2Ccenter 2C surrounding _ 2Ccenter (13)
= +
A’ DC 2 DG 2

We had to be careful in choosing the size of the grid spacing. The explicit
method goes numerically unstable if the grid spacing to time step ratio is too small.
Since we have already chosen a time step, we wanted to pick a grid Spacing that was
small enough but not be unstable. Substituting these expressions into equation 11 yields
a result independent of the concentrations at the center and the 4 closest points.

The ratio of grid space squared over time step has to be greater than 5 times the
diffusion coefficient. We determined that 24 micron grid spacing that we used was

numerically stable by using equation 15.

§D_N(4)gl <14)
4062
DC 2 JSDAI (15)

For a 3D grid the criteria for numerical stability is expressed by equation 16. The
coefficient under the square root is the number of nodal points surrounding one nodal
point plus the nodal point at the center (i.e. 4+1 for 2D and 6+1 for 3D). The critical grid
spacing for 3D diffusion is only 18% larger then for 2D diffusion. This expression will

be helpful when the model is expanded for diffusion in 3D.

00 2 77m: (16)

23

An alternative method of calculating the diffusion may be employed to maintain
stability at smaller grid spacing. A method for calculating the gradient a bacteria sense’s
independent of the grid spacing may be implemented based on the time it takes for a cell
to notice a concentration change. The time it takes a bacterium to react to a concentration
change is on the order of a run duration, which for E. coli is around 4 seconds (Berg,
1993).

3.4 ASSUMPTIONS

A cell swims in a series of runs and tumbles and randomly picks a new direction
based on the turn angle distribution. In this model cells are assumed to swim in a straight
line between tumbles and during a tumble only the direction of the cell changes and the
location remains the same.

Cell growth and death does not occur in the model. The number of cells remains
constant, which reasonably models cells in the lag phase (no growth and no death), and
the stationary phase (growth equals death). Migration for the exponential growth phase
may be modeled for simulations shorter than the cell doubling time.

All the cells are assumed to have the same swimming velocity, basal tumbling
probability, chemotactic sensitivity, and disassociation constant.

Multiple cells can occupy the same location on the grid and cells cannot sense the
presence of other cells directly. Only through the gradients in the chemoattractant
produced by the other cells can the neighboring cells be sensed.

Multiple gradients are not sensed, effects from other chemoattractants such as
oxygen are not currently modeled. The substrate used for energy in the nonconsumable

attractant simulations is assumed to be plentiful and not be a chemoattractant.

24

The availability of electron donors or acceptors and its effect on the consumption
rate is not considered but assumed to be available in sufficient quantity.

3.5 DIMENSIONALITY

The input parameters for chemotactic sensitivity and single cell swimming
velocity are three-dimensional values, which are scaled to two-dimensional values when
studying two-dimensional migration with consumption.

Equations 17 and 18 scale the dimensionality of the cell motion. The chemotactic
sensitivity is proportional to velocity squared. The differential tumbling frequency is
denoted as 1), and NT is the total number of receptors (Frymier et al., 1993; Rivero et al.,
1989). For simplicity the two-dimensional chemotactic sensitivity and velocity will be

noted as x0 and v respectively in this paper.

v21) : v30 Q (17)
J3

2
] (18)

2D_ 2 _ 3D

SIS

4. LINEAR GRADIENT NON-CONSUMABLE

The simulation parameters used in initial simulations were 2000 cells which is
small compared to the number encountered in experiments, but large enough to get
realistic average results without having to use excessive computer time (some 1 hour

simulations required 4 hours using an 800 MHZ Pentium).
Other authors have used a chemotactic sensitivity (3(0) value of 4.1E-4 cm2/s for

E. coli responding to or—methylaspartate (Strauss et al., 1995). A X0 value of 2E-2 C1Tl2/S

25

was selected for our simulations to give a noticeable chemotactic response with a small
linear gradient of 1.25E-5 mM/um that could be experimentally obtained easily. The
larger value of 10 is used solely to Show the trend in motility when cell motion parameters

are changed. Furthermore, variability in cell cultures can cause the X0 value to vary 20 %
from the mean (Marx and Aitken, 2000). The chemoattractant response is also heavily
dependent on temperature. Experiments by Adler showed that increasing temperature
from 20 C to 30 C yielded a 20—fold increase in bacteria attracted and no chemotaxis
below 15 C (Adler 1973). Since the X0 is dependent on the velocity squared, by only
doubling of cell velocity could explain a 4-fold increase in the chemotactic response.
Above 37 C the synthesis of the flagella stops which caused the movement of the bacteria
to diminish (Adler and Templeton, 1967).

4.1 SINGLE CELL TRAJECTORY

Figure 6 shows the simulation trajectory results of a cell movement according to
random motility compared to cell movements demonstrating chemotactic behavior using
a fixed, linear, non-consumable gradient, increasing in the positive y-direction. The
symbol chi represents the chemotactic sensitivity. The Simulations for the trajectory
results were completed using various chemotactic sensitivity coefficients. The velocity
and tumbling probability of the cells were fixed at values of 22 It m/s and 1.2 s",
respectively. Table 1 contains the input parameters used in the simulations of cells
responding to a linear non-consumable gradient.

When X0 is at or approximately zero, the cell clearly demonstrates a typical

random motion with very little migration away from the point of inoculation. While

gradually increasing the value of x0, the cell movement deviates from random motion to

26

that of a cell responding to a chemical gradient. If an unreasonably high value of X0 is
used, such as 4.1 cmzls, chemotaxis is clearly evident because the cell displays a very
strong, biased movement along the gradient without deviating from its path. The average

distance a cell moves from its initial location is larger for larger values of the chemotactic

 

 

 

 

 

 

 

 

 

 

sensitivity.
200 ‘
‘\ x0 (cm2/s)
A 0.02
E
a 50 «
>..
0 _. l __ _____F____a
-50 - .._._____..i_.._
-100
-200 -100 0 100 200

Xfum)
Figure 6. Cell trajectory comparison

Table 1. Values used in the simulations of a non-consumable linear gradient of a-methylaspartate for

E. coli.

 

 

 

 

 

 

 

 

 

 

NP 2000

x 40,000 um
Y 40,000 um
DG 40 um

p0 1.2 s‘I

v 22 urn/s
C 1 mM

K, 0.125 mM
103” 2E-2 cmZ/s

 

 

 

27

4.2 EFFECTS OF VARYING CHEMOTACTIC SENSITIVITY COEFFICIENT

Figure 7 shows the difference between cells with a moderate (Xo = 2E-2 cmzls)
chemotactic sensitivity and those demonstrating no chemotaxis at all (X0 = 0 cm2/s).
The lack of a gradient within the environment will cause the cells to move in a random
walk, thus denying any movement biased toward one direction. The cells will inhabit the
location close to the point of inoculation. An introduction of a gradient will cause the
cells to migrate away from their origin. Even at small values of X0. the cells will

demonstrate chemotaxis, biasing migration more with larger gradients and at longer time

periods.

 

5000

5

.3-3 . .L 3'1, ; 0.02(cm2/s)

2500 -

 

 

Y (pm)
o

 

-2500

 

 

 

-5000 r r i

-5000 -2500 0 2500 5000
X (um)

 

Figure 7. Snapshots of chemotactic and non-chemotactic cells in a linear gradient.

When the chemotactic sensitivity is small, they stay in a circular pattern and the
population moves to areas of higher concentration. When the chemotactic sensitivity gets

larger, they tumble less frequently and the shape of the colony flattens out to an oval

28

shape. At a very large chemotactic sensitivity a limiting value is reached where the cells
don’t tumble and they continue to move up gradient in the direction of their initial
orientation. The cell density profile also shows this result. Figure 8 is a comparison of
snapshots showing cell positioning at several values of the chemotactic sensitivity
coefficient (x0 = 4.1, 4.1E-1, 0) cm2/s taken at 30 minutes after inoculation. While
increasing X0. a cell population will migrate along the y-axis in the positive direction. At
low values of X0, there is little, noticeable movement of cells in the positive y-direction,
imitating non-chemotactic behavior at small time periods. If x0 becomes extremely large
(e. g. X0 = 4.1 cmzls), the cells cease to frequently tumble. AS a result, they move radially
in the direction of the gradient producing a cell snapshot resembling a semicircle, where
very little random motion occurs. The density of cells increases at locations of higher
chemoattractant concentration, proving that chemotaxis displaces cells more dramatically
than random motility.

To better illustrate the difference between the average displacement (drift
velocity) and the motility coefficient (diffusive mean square displacement) in Figure 8
the motility coefficient is larger for the cell distribution with a X0 of 4.1 cmz/s then the
cell distribution with a X0 of 0.41 cmZ/s, while the average displacement is larger for 10 of
0.41 cmZ/s then X0 of 4.1 cmZ/S. In different applications the goal may be to disperse the
cells more or it may be to transport the cell to a focused area so both the motility
coefficient and the average displacement of the cell population are important in defining

the behavior of the cell population.

29

 

50000

 

 

 

 

 

xo (cmz/S) 0
40000 ' 0'41
-4.1
30000 4"“

 

20000 // . . .\
10000 '37: ' ' "

 

Y (um)

 

 

 

 

 

 

-10000 I T I T r T F
40000 -30000 -20000 -10000 0 10000 20000 30000 40000

X (pm)
Figure 8. Snapshots of cell location along a linear non-consumable gradient.

Simulation results for the density profile of a 2000 cell population migrating in
response to a fixed linear chemical gradient are shown in Figure 9. The cell density is
plotted as a function of location in the y-direction for cells using a variety of X0 values.
Again, chemotaxis can be seen from observing the density of cells with respect to the
distance the cells traveled, as plotted on the x-axis of Figure 9. As the value of X0
increases from zero to 4.1 cmZ/s, the spatial location of the cells increases as well, clearly
showing the behavior of chemotaxis. Beyond a x0 value of 4.1 cmzls, the motility ceases

to increase, demonstrating a limit to how large the motility can get.

30

 

 

 

 

 

 

 

 

 

 

 

 

 

2 + 0
50 q j x, (cm ls) + 0.02
. ---— 0.2
3}; —~— 0.41
a 40 - '
3 ll
"6 30 -
l-u i . '
.3 i 1 l. .T
s 20 i T" '1
z i 1b 8
10 i ’ f
0
-5000 5000 15000 25000 35000 45000
Y (um)

Figure 9. Cell density profile in a linear gradient at varying chemotactic sensitivity coefficients

Figure 10 shows the effect of the chemotactic sensitivity coefficient (3(0) on the
motility and the mean displacement of the cells. The dotted line corresponds to the mean
distance up gradient at its respective 10 value and the solid line corresponds to the
motility data (displacement from its original position). The motility coefficient increases
as a sigmoid curve as plotted on a logarithmic scale, it starts to level out around X0 of IE-
] cm2/s for simulations using a fixed linear gradient (1 .25E-5 mM/um). The chemotactic
response increases exponentially over a range of X0 for approximately 4 orders of
magnitude. A dramatic increase in motility occurs at this critical sensitivity. A limiting
motility coefficient, corresponding to straight-run cell motion, is attained at very large x0
= 4.1 cmzls values. The calculated motility coefficients only correspond to the conditions
and simulation time used for determining them. With a different X0 or a different
concentration gradient, the motility will correspondingly change. The mean displacement

of the cell population increases with respect to an increasing chemotactic sensitivity

31

coefficient and levels off at higher chemotactic sensitivities are reached where cells in

 

 

 

 

 

 

 

 

 

 

 

frequently tumble.
1E+0 1 . , - 25000
20000
A 15000
Q
as 10000 a;
3 5000 g:
a 0 :1
a e
2 -5000 a
~10000
~15000
-20000
1E-8 r i r i T -25000

 

 

 

113—5 lE-4 113-3 1E—2 lE-l 1E+O 1E+1
Chemotactic Sensitivity Coefl’icient (cmz/s)

Figure 10. Chemotactic Sensitivity versus Motility and Mean Distance.

4.3 EFFECTS OF VARYING CELL SWIMMING VELOCITY
The cell swimming velocity parameter appears in equation 1 to suppress tumbling and
calculate the distance a cell moves to each new location. Thus, the effect that velocity
has on the motility coefficient is much more complex than the remaining parameters.
According to equation 18, the chemotactic sensitivity is a function of velocity so it needs
to change as the velocity changes and should not be kept constant. To determine the
effect velocity has on motility the DNT was kept constant. A typical value of UNTIS 75 5,
however a value of 750 s was used in the simulation to obtain a more pronounced
chemotactic effect (Ford, 1992). Simulations were completed at several cell velocity
values for both the random motility and chemotactic scenarios, in the presence of a fixed,
linear, methylaspartate gradient. Figure 11 gives the trend of the motility coefficient

calculated for a few of the velocities.

32

 

 

 

 

 

 

 

 

 

ASE-5 +motility 0 i DNT = (S) o’ 4500
4013-5 — - -o - motility 750 - -—'- * “ __- .‘T‘ "i 3000
3.535 4 +position0 l i- .7 ‘ J—:«— +—--
- -o - position 750 j . . ' .’
3 3.05.5 -- -_-___-_ -- — -__.__ '—:—-’T'“"‘ ' “i 1500
E
3 0
5
e -1500
E
-3000
—4500

 

 

 

 

Cell swinurl'ng velocity (um/s)

Figure 11. Velocity versus Motility and Mean Distance.

The velocity ranged from 1 uni/S (chemotactic motility 2.65E-9 cmzls) to a much
higher velocity of 70 urn/S (chemotactic motility of 4.2E-5 cmzls). The circles represent
mean distance results whereas diamonds represent the motility results, while dotted lines
refer to chemotactic data and the solid lines refer to non-chemotactic.

The dotted line with diamonds displays the motility trend with respect to
increasing velocity at a constant UNT of 750 S. The motility increases in a quadratic order
for both the random and chemotactic cases, as predicted by equation 2 for random
motility. The mean displacement of the cell population from their initial location
remained constant for non—chemotactic cells as velocity increased. For the chemotactic
scenario, the mean cell displacement increases as a quadratic, analogous to the non-

chemotactic case.

33

4.4 EFFECTS OF VARYING THE TUMBLING PROBABILITY

Figure 12 demonstrates the effect of adjusting the tumbling probability of the cells
for both the chemotactic and non-chemotactic scenarios. The chemotactic sensitivity
coefficient X0 of 2E-2 cmzls was used with a constant cell swimming velocity of 22 urn/s.
At 30 minutes of cell simulation, the center of mass displacement remains constant as the
tumbling probability is increased for both chemotactic and non-chemotactic cells. The
simulation results Show that chemotactic and non-chemotactic populations will disperse
more at lower tumbling probabilities. The chemotactic cells will be displaced up gradient,
while the non-chemotactic cells will remain clustered around the origin. Increasing the
tumbling probability will cause the cells to tumble frequently, resulting in little dispersal
of the cell colony. For both chemotactic and non-chemotactic scenarios, the motility
decreases dramatically until it reaches a tumbling probability of 1 s", which is common

for E. coli, and then it levels off.

 

 

 

 

 

 

 

 

 

 

6E-5 _ _ _ _ _. _________ , ________ 3 2500
2 —~ 2000
.. rain.
I» - - - 0t - — 1000 A
..\ 413-5 ~ +Position0 E
5, - -o - Position .02 - 500 I?
.3‘ O -
g -500 E
2 -1000 E
-1500
—2000
-2500
0 1 2 3 4 5
mums probability (s")

Figure 12. Tumbling Probability versus Motility and Mean Distance.

34

For Simulations running at very Short time intervals (e.g. a few seconds or less),
authors have reported that cells with shorter run lengths (higher tumbling probabilities)
Show more chemotactic response, thus contradicting our simulations results (Dillon et al.,
1995). We have found that chemotactic response (mean displacement) at longer time
intervals (e.g. 30 min) was independent of the tumbling probability. It will take longer
for a cell with extended run lengths to orient itself with the gradient, Showing a
chemotactic response from its initial location. A random orientation of the cell, if facing
down gradient, will move with random motility until it readjusts itself in the direction of
the gradient, thus allowing it to bias its motion. A longer amount of time is required for a
cell with extended run lengths to complete enough tumbles to orient itself with the
gradient.

At shorter time intervals, cells with high tumbling probability will often orient
themselves more rapidly, allowing movement up gradient. A shorter length of time is
then required for movement up gradient, with the cells oriented in a tight cluster. At
longer run times, chemotactic response (mean displacement) of high and low tumble
probability is indistinguishable, but our results conclude that cell motility decreases with
higher tumbling probability due to the reduced dispersion of the cells (Figure 12).
Simulation results Show that the motility coefficient increases at low tumbling
probability.

4.5 GAUSSIAN TURN ANGLE DISTRIBUTIONS

A set of gaussian turn angle distributions were calculated and normalized to

determine the effect the mean turn angle and the variance of turn angle distribution would

have on the motility coefficient. These distributions are symmetrical with both the mean

35

and the median at the peak of the distribution. The gaussian turn angle distribution was

calculated using equation 19.

_ 2
9Mean )
2
‘7 (19)

 

probability = Exp —

Figure 13 and Figure 14 Show the normalized turn angle distributions used in our
Simulations for the mean value and the variance, respectively. The turn angle distribution
used for all the other Simulations in this thesis is for E. coli strain NRSO reported by

Frymier (Duffy and Ford, 1997).

 

 

 

 

 

 

 

 

 

 

0 3O 60 90 120 150 180
Mean angle (degrees)

Figure 13. Normalized turn angle distribution varying mean.

36

 

 

 

 

 

 

 

 

 

 

0 30 6O 90 120 150 180
Varience (degrees)

Figure 14. Normalized turn angle distribution varying variance.

4.6 COMPARISON OF E. COLI AND P. PUT IDA MOTILITY

The actual turn angle distributions are not gaussian but are skewed toward smaller
turn angles as shown in Figure 15. The gaussian turn angle distributions are symmetric
thus they have the same mean and median values that occur at the peak to the
distribution. The turn angle distributions for E. coli have mean, median and peak values
distinctive from each other, while P. putida has a bimodal distribution with two peaks.
When correlating the behavior of the bacteria with different turn angle distributions,
which value contributes most to the motility behavior of the bacteria? One method would
be to look at the median turn angle, but the correlation would not be very smooth for
skewed distributions. A better method would be to take the moments of the distribution,
which have a strong Statistical basis. The first moment would be the mean or average of

the turn angles, while the second moment would be the variance. According to our

37

simulations, the mean of the turn angle is the most important mode in determining the
motility coefficient. Table 2 shows the values of the mean, median and peak of the cell

turn angle distributions.

   
   

 

- - - - P. putida
_ ____ E. coli NR50
—6— E. coli AW405

-—.- ‘777___a—_-_7_.--_4

 

 

 

 

 

 

O 20 4O 6O 80 100 120 140 160 180
Tumangle (degrees)

Table 2. Mean, median, and peak values of bacteria turn angle distributions.

 

 

 

 

 

Cell type Mean Median Peak
E. coli strain NR50 79.2 75.7 68

E. coli strain AW405 65 59.2 50

P. putida 85 77.2 35, 145

 

 

 

 

 

Figure 16 Shows a comparison snapshot of the cell density of E. coli and P.
putida. The figure clearly Shows that differences in the shape of the turn angle
distribution do not Si gnificantly affect the motility. Even though our simulation method
was developed for peritrichious bacteria similar to that of E. coli and Salmonella, it works
equally as well for Pseudomonas putida. Pseudomonas putida has polar flagella, a
bimodal turn angle distribution and swims twice as fast as E. coli (Duffy and Ford, 1997).

Even though P. putida swims twice as fast as E. coli, all the parameters in the Simulation

38

snap Shot in Figure 16 where held constant except the turn angle distribution for

comparison reasons. Comparisons of the turn angle distributions are in Figure 15.

Y (um)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

7000
° P. putida
6000 .E. coli
5000 tiff... if:
s I" t ”-u is: e; "alt-k:
. O‘d' ...v‘d' or“?
4000 "51% ”CL“?! .'.' ..
I! we ‘ 51"“?ch 6 11.":
3000 e .‘r :‘Efififi; 53":
J . P} . ."
20““ T's: W?» :W
1000 T , ‘_1'5,«'tj -’
0
-1000 , 1 '
-4000 -2000 0 2000 4000
X (pm)

Figure 16. Cell density profile snapshot of E. coli and P. putida at 30 minutes.

4.7 EFFECTS OF THE MEAN AND VARIANCE OF THE TURN ANGLE DISTRIBUTION

The variance of the gaussian distribution was fixed to a value that mimics that of

E. coli (30 degrees) and the turn angle was adjusted to values of 0, 45, 90, 135 and 180

degrees. A value of 0 degrees corresponds to the cell continuing to move in its present

direction whereas 180 degrees refer to cells making a full reversal. While using a x0

value of 2E-2 cmz/S, motility and mean displacement data was taken at these various turn

angles. Figure 17 clearly shows that cells utilizing smaller turn angles will migrate

farther than those with turn angles that are large. Thus, smaller turn angles are desirable.

39

 

2.5E—5 ‘ _____ _' _____ ,‘ _____ 5 _____ '3 2500

 

 

 

 

 

 

 

\ x0 (cmzls) T 2000

2013-5 ~ “\ Tflmfly}, 002— 1500
"‘ \ ' ' 0t - 1000 A
cg ‘\ +Position0 E
g 1515-5 4 . - a - Position of 500 5
5‘ . . o 0 g
tawny ’ """ '9""' -500 '5
e e
z -1000 9‘

5.0E—6 ~ -1500

-2000

0.0E+O T . T I I 1 -2500

O 30 6O 90 120 150 180
Mean turnangle (degrees)

Figure 17. Turn angle mean value versus motility and cell mean distance in a linear nonconsumable
gradient.

Figure 18 shows the effect of the turn angle distribution variance on both the
motility and the mean displacement of the cells. Simulations were run for both the
chemotactic and non-chemotactic scenarios. The variance was adjusted while
maintaining a turn angle of 90 degrees (right-angle turns). Using a X0 value of 2.0E-2
cmz/s for the chemotactic case, the turn angle distribution variance was adjusted to 10,
30, 60 and 180 degrees. The process was then repeated for the non-chemotactic scenario.
Adjustments in the turn angle distribution variance had no effect on both the motility and
the mean displacement of the cells for both scenarios. That is, the motility remained

constant at every angle tested as well as the mean displacement of the cells.

40

 

1.235 .__3_._. _______________ $2500

9 - 4- - . .............. .
LOB-5 - ’ x0 (mo/S) .

+ Motility o

 

 

 

 

 

 

NE? 8.0E—6 ~ - 4 - 11:40tility.(())2 ” 1000 E
+ osrtion _
E - 4 - Position .02 500 5
3. 6.0E—6 . H—+ 0 =
O T” -500 8
2 4.0E—6 - __ _1000 9..
2.0E-6 ~ ‘" "500
._.__. T —2000
0.0E+O . r t . , t -2500

 

 

0 30 60 90 120 150 180
Variance (degrees)

Figure 18. Variance of the Turn Angle versus the motility and cell mean distance in a linear non-

consumable gradient.

The accuracy of our model was checked against the widely accepted theoretical
equation for random motility (equation 2). This equation accounts for tumbles as
instantaneous so simulations with bacteria swimming runs at each time step were
performed to make a true comparison. The results of our simulation show excellent
agreement in motility for both 2D and 3D cell movement with the motility calculated
with the theoretical equation. Plots of motility vs. mean turn angle (Figure 19) and

variance (Figure 20) show good agreement over the entire range of angles.

41

 

 

—-A— 2d, 3d eqn
" ‘ - o - 3d sim
A" —0— 2d sirn

 

 

 

7_—7.~v—-<t

 

 

 

 

 

 

O 20 4O 60 80 100 120 140 160 180
Mean turn angle (degrees)

Figure 19. Comparison of the mean turn angle in 2D, 3D against theoretical equation.

 

 

 

 

 

 

 

 

4.0E—6
+2d 3d eq

3.5E-6 n.-- ~ *7 ~ — ~ ~ ~ ~ —-~~--~— — —- w—g— . —0— 2d sim
7; 3.0E-6 arm» —-- _ _,_, ~ ~ I- ”I 4 ,-_-_ _ _.. _ 3d sim
"5 2513-6 we , —* -_ ~* ~ # --___-___,,_.__ —r— — _ _ ._ ~ _- _ .. -...,_ __ _ __
ii
5. 2.0E—6 4 - ~ ~~ , , , ~ ‘ - , — ~——~—-—-—»-_.—_—-___-_---_“CACHE -___.
so. _ ____
2 LOB-6 — ‘_~ _,_ —

5.0E-7 -»— , , , I - - — -.I-__I — ——--————_-_.-__-mm-.- ,___-_.I.~___

0.0E+O t t r r , , , I

 

0 2O 4O 6O 80 100 120 140 160 180

Variance (degrees)

Figure 20. Comparison of the turn angle variance in 2D, 3D against theoretical equation.

42

5. CONSUMABLE NON POROUS

Knowledge about chemotaxis in a fixed linear gradient can be used to better
understand chemotaxis in a dynamic chemoattractant concentration environment. The
environment within the presence of a chemoattractant is constantly being depleted based
on the location of the cells and then simultaneously replenished by diffusion. This results
in a complex concentration profile, which can be approximated as a fixed linear gradient
at each cell location for the duration of one time step, thus resulting in a variety of
gradient directions. The consumption coefficient (k) is a measure of the rate of a cell’s
substance intake through its cell wall, which is carried out by transport proteins located in
the cell wall. The rate of uptake is affected by the saturation of the receptors by the
concentration around the cell (Madigan et al., 2000). The rate is also affected by the
availability of the electron donors, acceptors and the cell metabolism (Dillon and Fauci,
2000).

Consumption of chemoattractants causes a band or wave of cells to migrate from
the center. Multiple bands of cells have been observed from experiments (Figure 21) and
are most likely the result of consumption of multiple chemoattractants (Emerson et al.,
1994). Table 3 contains the parameters used in the simulations where cell react to

gradients created from consumption of the chemoattractant.

43

 

Figure 21. Multiple Chemotactic waves of bacteria due to consumable chemoattractant.

Table 3. Values used in the simulations starting with a uniform concentration of aspartate consumed

by E. coli.

 

44

Figure 22 demonstrates how cells migrate outward from the point of inoculation

in response to local attractant gradients formed by consumption with a chemotactic
sensitivity of 113-5 cmzls. In the snap shot of consumption rate of 0.1 s'lcell'l the

gradients are too small to bias all the cells away from the center. At the consumption rate
of 1.5 s'lcell‘l all the cells will detect a large gradient and will bias their migration away
from the center. At the larger consumption rates (2.0, 3.5) s'lcell'l, the gradients are
depleted at the center and a portion of the cells will move with random motility while the

cells at the perimeter will continue to bias there migration away from the center.

 

Figure 22. Snapshots with Consumption showing the cell density patterns.

The value of the consumption rate that each pattern happens is different for each
set of experimental parameters. The number of cells, initial concentration of

chemoattractant, the swimming parameters of the bacteria and the chemotactic sensitivity

45

and the consumption rate will effect when each pattern will be produced and do not
necessarily correspond to values used in other parts of the paper.

The radial cell density shows that at higher consumption rates, cells can be
stranded in areas of depleted chemoattractant and will move with random motion while
the cells on the outer band will continue to consume the chemoattractant and migrate
further away. Leaving cells behind lowers the motility of the population. This also
explains the bands formed during experiments. Only 3 of the 4 cell density profiles are

shown to enhance clarity. The 30~minute simulation results in the Figure 23 used a

chemotactic sensitivity of lE-3 cmz/s for consumption of aspartate.

    
   

5.0E—4

 

_ [ Dimensionless concentration )fi _ 0-1 -"- 1-5

 

 

 

Number of cells ~ sec onds

 

 

 

 

 

 

 

 

 

 

0 500 1000 1500 2000 2500 3000 3500 4000
Radial Distance (pm)

Figure 23. Radial cell density at different consumption rates.

5.1 INITIAL DISTRIBUTION

Simulations were completed with cells initially located at every grid point, thus
forming an equal distribution of cells. Within the presence of a constant concentration of

chemoattractant, an increase in the consumption coefficient or x0 does not affect the

46

motility. The concentration is constantly being depleted at every point in equal amounts.
Therefore, no significant gradients will form. The random movement of the cells will
cause small variations in cell density. As a result, small variations in chemoattractant
concentration will occur with a random direction Of the gradients.

5.2 DIFFUSION

Diffusion occurs naturally from a difference in concentration (concentration
gradient) between two nodal points. The rate of diffusion is proportional to the gradient
size. A linear gradient is achieved when the diffusion into and out of a point is equal.
Inoculation of all the cells at one location produces a cooperative gradient in the radial
direction from that location, giving different results than an equal distribution of cells.
Since the gradient is radially symmetric, the cells will produce a gradient in the initial
direction the cells are facing.

Figure 24 Shows that there is an optimum consumption rate that maximizes
motility. This corresponds to the point where the produced gradient is the largest without
leaving cells behind. In Figure 24 the motility increases to a peak value, which occurs
from balancing large gradients with maintaining gradients at the center. At low (0.1 s’l
cell") and high (8 s'lcell'l) consumption coefficient values, the motility approaches that
of random motion, increasing to a peak value, but has a slower drop off. The random
motility coefficient calculated in the absence of a chemoattractant gradient using a 0.1 s
time step is 1.32E—6 cmz/s. The magnitude of the peak is proportional to the x0 value.
Widman also observed that the carrying capacity of the gradient with consumption is

larger at increased chemotactic sensitivity values (1997).

47

 

 

ZOE-4

2 +5E—6
xo(Cm /S) +SE-4
. -4 - 7 -7. 7 7 -- g — —————__..— -)l(-1E-3
1 6E +5E-3

 

 

 

 

 

 

 

 

 

O l 2 3 4 5 6 7 8
Consunption Coefficient (S'lcells'l)
Figure 24. Effect of the consumption coefficient on motility with diffusion.

At consumption coefficients less than 2 S'lcell", smaller gradients exist, but the
chemoattractant is not completely depleted because diffusion occurs equally as fast as
consumption. The cells subsequently can find an existing gradient and begin migration
from the center of inoculation. At consumption coefficients greater than 4 s"cell'l, the
existing gradients are larger, but the chemoattractant is completely depleted at the center,
thus a majority of the cells are left behind and move with random motility. This causes a
decrease in pOpulation motility. The consumption simulations correlate with the error for
the linear gradient also yielding less than 3% error. The error bars are smaller than the
graph symbols, so they were left off for clarity.

5.3 NO DIFFUSION

The most noticeable difference between the diffusion and non-diffusion

simulations is the diffusion case has a higher peak value (Figure 25). At X0 of 5E-3 cm2/S

48

the maximum motility coefficient calculated to be 1.34E-4 cmZ/s for diffusion and 2.4E-5
cm2/s for the no diffusion simulation. The diffusion simulations have dark lines and no
diffusion has gray lines.

The results Of these simulations indicate that diffusion is important in cell
migration with consumption and Should not be ignored. For non-diffusion simulations,
only the cells at the front of the chemotactic wave will experience gradients. These cells
will depleted the chemoattractant completely leaving most of the cells at the center
moving with random motion. Increasing the chemotactic sensitivity increases the

motility for all the values of consumption coefficient that were simulated.

 

 

 

 

 

7'05'5 x (mg/S) +554 No diffusion
60135“ f _ ° g +553Nodifl‘usion
- ' “FTC F “E _ +5E—4Di1fuson
+5E-3Dimlsion

E mm5~~——m --n «—~w~—~~“-——»

E Mm6~~+-~~e --’_~~ -4—4——~

.39

g 3.055 —t 4» ~ ~—-~~ —«

G . .

2 2.055— .~—~——u - —--—————————--~——~—

 

1.055 ,7, — — _. A- _ .__._-__ -_._.__-_ ____.-

 

 

0.05+0 ' t . t . t , .

O 1 2 3 4 5 6 7 8
Comunption coefficient (S'lcells'l)

Figure 25. Motility with diffusion and no diffusion on a logarithmic scale.

An interesting phenomenon is at low consumption coefficients the non-diffusion
case has a higher motility coefficient because when none of the cells get left behind in
either case the only effect of diffusion is smoothing out the gradients. At higher
consumption coefficients the diffusion case has a higher motility coefficient because the

gradients are maintained and the cells don’t get stranded. The effect of diffusion on the

49

motility coefficient for the equal distribution of cells doesn’t have much of an effect
because of the lack of gradients to smooth out.

5.4 ERROR ANALYSIS

Motility values were obtained by running replicate simulations, while altering the
seed number for the random number subroutine and keeping all other parameters
constant. A standard deviation was calculated for each scenario and used to determine an
error percentage. The percent error for all simulation results, including the linear
gradient and consumption cases, is at or below 3%. A notable trend with regard to the
error analysis is the percent error increases with an increase in motility. At larger
motility values, the cells will move farther from their initial positions, giving them more
possible locations thus more variability.

5.5 ALTERNATIVE METHOD OF GRADIENT DETECTION

The decaying step function is ideal for validating chemotactic models with
consumption of a chemoattractant. The decaying step function creates both changing
temporal and spatial gradients and is easily defined. Comparisons of different parameters
can be done knowing that the gradient will be the same for each case. The concentration
profiles in each case are similar resulting in a common gradient detection response.

The simulations start with a uniform concentration of 0.1 mM aspartate and all the
cells starting at the center yielding a plot of the cell density and the concentration profile
due to consumption of the chemoattractant. The peak of the cell density is located at the
lower concentrations of the formed gradient. In Figure 26, the chemoattractant
concentration is completely diminished resulting in a second peak of cell concentration

located at the origin that moves with random motility. In Figure 27 almost all the cells

50

are biased away from the center resulting in a wide peak compared to plot a where a

narrow peak is separated from the center peak by an area of low cell concentration.

 

 

 

 

 

 

 

 

2.55+4 0.125
S

e 2::
¥ o
@1554 0.075 g
ii

1.0E+4 0.05 o V
"3 5
= 5
55.059‘ 0-01’5 g

0.05+0 - -- . 0

 

 

O 20 40 6O 80 100 120
Radial distance from innoculation (um)

Figure 26. Cell density profile where chemoattractant is depleted at center.

 

 

 

 

 

 

 

 

 

 

 

2.55+4 , 0.125
r.-
1 -.-.-
A2.05+4 0,1 g
‘ 8
€1,554 0.075 a A
3 ..
a? > 5 E
51.0154 ; 0.05 g
g <— g
05.059 . 0.025 2
\yu U
0.0E+0 I I T 1‘ I O
0 20 40 60 30 100 120

Radial distance from innoculation (u m)

Figure 27. Cell density profile where chemoattractant is not depleted at center.

51

An alternate method of detecting the chemoattractant gradient has been developed
to account for temporal gradients and compare to Spatial gradients, since experiments
have indicated that bacteria react to changes in concentration over the 4 previous seconds
(Berg, 1993). In the case of the decaying step function, the center position of the
concentration gradient is fixed and the cell density changes as it reacts to the gradient.
The gradient also decreases with time. Hypothetically cells at the top of the gradient will
detect a lower concentration temporally and may chemotactically go the wrong way
down the gradient. Cells at the low end of the gradient will sense a higher concentration
over time and may swim away from the gradient. This will not be a problem because
simulations Show cell density bias up gradient within 1 second. A comparison of
different chemotactic sensitivities on the cell density in response to a fuctose gradient at 6

minutes is shown Figure 28. Initially the relative cell density is uniform at 1 all along the

 

 

 

 

chamber.
1 0.25
-)(— 1.05E-02 x0 ’(cmzls)
— 3.50505 '" ’ 0 2
—Chemoattractant- '

 

 

r

l
.9
L7.

 

 

Chemoattractant concentration
(mM)

 

Relative cell density
0 -— N t»
O M '— UI N U! U) Ut &
L l I I L l
I
l
l
l
1
1
l
l

 

 

-4000 0 4000
Chandler position (pm)

Figure 28. Comparison of different chemotactic sensitivities on the cell density profile.

52

Simulations using the values for fuctose were used to compare our simulation
results to those obtained by Ford’s research lab for cells responding to a decaying step
gradient. The shape of our relative cell density profiles correlate well with those obtained
by Ford’s research group.

The original mathematical relation for the chemotactic biasing of the tumbling
probability included a time derivative as well as the spatial dot product. Ford assumed
that the time derivative was much smaller than the spatial derivative so it was dropped

from the simplified equation.

30
lo Kd (161;?

v (Kd+C)2 v dt

 

pzpoexp— +SOVC)
(20)

Our simulations show that the time derivative is on the same order of magnitude
as the spatial derivative for the simulations with a chemotactic sensitivity for fuctose of
3.5E-5 cm2/s. The combined simulations with the spatial and time derivative resulted in a
much enhanced cell density up gradient. Simulation using the full time derivative and
dot product results, indicated by the line with x’s, in a cell density profile much enhanced
over either one separately.

An alternate method of detecting the gradient by the concentration change at each

cell at 4 seconds as shown in equation 21 was developed.

10 Kd dC

V (I<d+6‘)2 [\[dX2+dY2 (2”

( K )dt24)

 

 

p=p0 6X1) -

 

 

 

 

53

The Simulation results in Figure 29 show a significantly lower motility for this
method where the data series with the squares represents the temporal gradient converted
to a spatial gradient.

Obviously this method under predicts the gradients that are seen using the regular
spatial gradient depicted by the regular solid line. The simulations using only the time
derivative result in a cell density profile depicted by circles that is slightly reduced from
the regular spatial derivative. Fuctose with a chemotactic sensitivity of 3513-5 cm2/s was

used in the comparison simulations at 12 minutes.

 

 

 

    

x, = 3.555 (cmzls)

spatial
—0- time derivative only
- ->( - spatial plus time drrevative
—13 —converted
0 t
-0.4 0 0.4

Chamber position (cm)

Relative cell denstiy

 

 

 

 

 

 

 

 

Figure 29. Relative cell density profiles of gradient sensing methods with x. of 3.5E-5.

A close correlation of cell density profiles resulted from the simulations using
105E-4 cm2/s for a chemotactic sensitivity for fuctose at 6 minutes for the entire set of
different gradient sensing methods except the simple temporal method, which didn’t

cause much biased migration as seen in Figure 30. The addition of the time derivative

54

does enhance the chemotactic migratiOn of the bacteria at lower chemotactic sensitivities,
but when the chemotactic sensitivity is very high, tumbling is already suppressed so the

addition of the time derivative will not have a significant effect as seen in Figure 30.

 

 
 
  
    
  

 

 

 

 

 

 

 

 

 

 

4
—regular spatial
3-5 A +tirre derivative only ——" — — — —«
3 _ — spatial plus tine derivative _
a —G— converted spathl
.3 2.5 .__~ ‘-—t~— ———~-~-
= x0: 10554 (cmzls)
3 2 EM. _ ,2- -2 L L .2 _ __2_- ..
3
'1: 1.5 -
.5
O
a: l ‘
0.5 4
0 l
-4000 0 400°
Chandrer position (um)

Figure 30. Relative cell density profiles of gradient sensing methods with x. of 105E-4
5.6 DEPENDENCE OF CONCENTRATION ON CHEMOTAXIS

The consumable chemoattractant concentration profile resembles a decaying step
function, but is shifted along with the bacteria as they migrate away from the center. The
peak of the bacteria density in Figure 31 stays at the base of the gradient resulting in the
continued sensing of a large gradient. The relation between the concentration profile and
the cell density profile that we obtained from the CD model agree with the experimental
profiles researchers have observed. Amperometric sensors have been used to measure
the spatial concentration gradients of the chemoattractants glucose and oxygen, while

simultaneously measuring the cell density with a CCD camera. In experiments with a

55

diffusion gradient chamber, the front of the chemotactic wave showed a high density of
cells at the base of the gradients of glucose and oxygen that were produced from

consumption (Peteu et al., 1998). In the simulation results without growth, the gradient
decreases as the bacteria thin out from covering a larger area, and diffusion replenishes

the area behind the wave front.

 

 

 

 

 

 

 

 

 

 

 

1E+6 . . . . 0.15
cell dens1ty 1.5 min —l-—cell densrty 6 rmn

9E+5 ~ - .0- - concentration 1.5 min - -I- - concentration6 min _ 0.13
"A 8E+5 “h 5* 5* ._-__ 011 A
E7E+5~k*—~— 'EE
‘3 —- 009 v
23 6E+5 -——*v~—~ - — — g ' E g
ESE-+5 l’ ”W 5* — ~— {~—~~-~—mfl-- ._- ~ «F 0.07 g E
'- 17___¥ i ,‘____.. I _ “-‘w-~ , _* ___ , a
= 3E+5 - , - “—JL— ~ ~‘—’ --— M“ ~
53 " ~— 0.03 ‘3

2E+5 -— —- —- ~— -—-~—— .

15+5 !~ 7 A ~- 0.01

OE ’3‘:5‘W1 “alums... -0.01

O 500 1000 1500 2000 2500 3000
Radial distance (pm)

Figure 31. Cell density profile and concentration profile due to consumption.

The concentration that the cell senses has an effect on how strongly the cell
population migration is biased when exposed to a gradient. A decaying step gradient of
methylaspartate was simulated with 3000 cells starting at the midpoint of the step
gradient. To determine the effect of the concentration, the concentration that the cell
senses was varied but the magnitude of the methylaspartate gradient was kept constant by
shifting the concentration up and starting with an initial step change in concentration of

0.2 mM. In Figure 32, cell density I responded to a concentration varying between 0 and

56

.2 mM, which shifted the peak up gradient 850 microns. Cell density 2 corresponds with
cells that had a step gradient of between .2 and .4 mM and shifted the peak 400 microns
up gradient. The peak of cell density profile 3 was Still centered at the cells initial

location responding to a gradient between .8 and 1 mM indicating saturated receptors.

 

 

 

250 , , +~ 4 — M ‘W‘T 1.2
—t—0-.2IIM
—.2—.4mM ’
200_ -—.8-lmM .2-__._ 1

 

 

 

 

 

 

‘1‘ w— 0.8
150 - -

 

 

 

 

 

 

 

 

E
a E
o I
s ‘* ‘ ‘***‘ i
E ‘ “F l 8
(I , g
50 4' —‘* ._:—7 0.2 a
O ......... M p O
-3000 -2000 -1000 O 1000 2000 3000
Positionqim)

Figure 32. The eflect of concentration on cell density in a decaying step gradient.

6. MIGRATION IN POROUS MEDIA

Migration of bacteria in pores is a widespread phenomenon and needs to be
understood to apply them to applications in bioremediation in soil and infection of
implants in the body. Porous media changes the dynamics of bacterial migration,
hydrodynamic effects of solid surfaces and sorption to the solid surfaces needs to be
taken into account. Figure 33 shows a schematic of the swimming path of a cell through

a bed of spheres. Just as important and not as widely considered is the production and the

57

magnitude of chemoattractant gradients can influence the direction and rate of dispersal

of bacteria.

[I
’,/” I .

I '4
I

I ’ r I--.:-. ’1
In

‘0 I,

I . ,

ll - 2!. a»: .
/,I,.v//,,,...’,,¢,.
“$9”??? "I; "r”

” i “a

Figure 33. Schemtic of cell trajectory through a pen or spheres.

 

 

Experiments with Pseudomonas KC reacting to a nitrate gradient showed farther
migration when swimming through porous soil core samples then when swimming in
bulk phase fluid (Figure 34). This result seemed surprising since one would speculate
that the migration in pores would be hindered due to collision and possible adsorption

with soil particles. The penetration ratio in this paper is defined as b/a.

 

Figure 34. Photograph of Pseudomonas KC migrating faster through sand core sample.

58

6.1 Ser PORES

The first set of simulations we did was look at the behavior of a population of
cells in a simple slit pore with and without consumption. In Figure 35, we set the initial
location of 2000 cells into the center of the slit pore of consumable chemoattractant.
After Simulations of 15 minutes, the cells indicated by the little black dots bias their
migration to the two ends of the pore. The chemoattractant at the center of the pore is
depleted indicated by the light gray color and a gradient points in the direction of the
openings of the pore as the cells follow the gradient. The black area at the Openings of
the pores indicates the high concentration of chemoattractant while the black diamond

indicates the average location of the cells.

 

Figure 35. Simulation results of cells in a slit pore.

We found that with a linear non-consumable chemoattractant, the motility did not

change very much as the pores got smaller. With the consumable chemoattractant, the

59

motility of the cells was more complicated. Under some conditions, the motility was

actually lower than random motility where gradients are not sustained and collisions with

the pore wall hindered movement as seen in Figure 36.

 

 

 

 

 

 

 

 

 

 

1.2E-4
_ Dimensionless concentration . 0'01
LOB-4 ‘ ‘ " "m " a Numberof cells-seconds ‘ 0'1
u -B—O.2

7; 8013-5 .. +fi A '1 a -~ - —~-~ ~ *7 ,2 — ,__,_. —)(—0.025
"E o x, = 5E-4 (cmzls)
3 6.0E-5 to " +~e ** ,_ ~ ~ , ~ ~ — - *~—--——-~- -5
3’
E40355 , ' *- —~— ~——-— if — -—-—
2 i, c

2055 ,3, ~— - .

', 1 \ :
00E+O —‘ t t t t t t
0 500 1000 1500 2000 2500 3000 3500

Pore width (pm)
Figure 36. Effect of pore width on motility at different consumption rates.

Figure 37 shows that there is an optimum population consumption rate for each
size pore. Cells in smaller pore sizes benefit from a smaller consumption rate. Larger
gradients could be formed quicker as the cells are confined until group consumption

outpaces diffusion so the chemoattractant is depleted and no gradient is sustained for

other cells to detect.

60

 

 

3.05—5 .
x0 = 5134 (cmZ/S) +50 I'B Wdth

—1 pore width
2513-5 ~ ‘-—---‘*—--~ r“ * +200 pore width
+400 pore width

1.55-5-+ i . ._5_ _ . - _ _
1.05.5 ~ , h— . a _ z-f-_4__s

V
5.0E-6 a n]. “a"- _ - -__ _ ‘i_

 

 

 

 

g
l
l

 

Motility (cmzls)

0
U

 

 

 

 

0.0E+0 ‘ . . 1 , .
0 200 400 600 800 1000 1200
Population cosunpfion rate (s'l)

Figure 37. Motility at different population consumption rates.

The combined term of the consumption rate times the number of cells gives the
population consumption rate which determines the chemoattractant gradient and thus
determines the motility. Different consumption rates with different number of cells but
the same product results in virtually the same motility indicating that a larger population
of cells could be effectively modeled by using a larger consumption rate and a smaller
number of cells. The motility for each population consumption rate rises to a peak almost
linearly and then drops off with a quadratic order at larger pore sizes (Figure 37). The
linear regime represents a sustained chemoattractant gradient and over time more and
more cells move away from the center and begin to bias their motion up gradient.

The curvilinear regime represents the continued depletion of the chemoattractant
at the center and only a fixed portion of the cells move chemotactically up the gradient

and the rest move with random motion. Smaller population consumption rates had

61

correspondingly smaller motilities with peaks at smaller pore widths and a smaller range
of chemotactic acceleration in pores. Larger cell populations have the ability to
accelerate their motility in larger pores giving them a competitive advantage against
smaller colonies of cells in the same size pore as shown in Figure 38. Smaller colonies
only have a marginal advantage in small pores since the gradient is depleted and only a

fixed amount of cells will sense the gradient.

 

 

 

 

3.0E—5
k = 0.1 Dimensionless concentration BE 3% 22g:
2-55-5 ‘ 0 Numberof cells - seconds +6000 C8115

  
  

x0 = 5E-4 (cmzls)

\

 

2.0E—5 5

 

Motility (cmzls)
27.
f.“
U!

l
P
l
l
l
l
l

 

 

 

Figure 38. Effect of pore width on motility with different number of cells.

To determine what causes enhanced migration in pores, a decaying step gradient
was simulated in different size slit pores (Figure 39). Three thousand cells responding to
gradients in aspartate with a chemotactic sensitivity of 5E-4 cmZ/s at 6 minutes were used
for these comparisons and a consumption rate of 0.1 s'1 cells'l was used for the
simulations of cells creating their own gradient by consumption. Smaller pore sizes
showed a small decrease in the motility for the random, and non-consumable step

gradient. The presence of a consumable chemoattractant enables accelerated migration in

62

pores. The magnitude of the acceleration is constant for larger pore sizes and increases

rapidly at moderate pore widths and then drops off at a critical pore width.

 

    

 

 

 

 

 

 

 

2.5E—S .
—9—non-consumable step gradient
2.0E—5 ~ *——--— __ _ " ' " ' “0 attractant
A —S—consmmble attractant
&
"51.513.5— _ ___-_.
8
.5‘
§ 1.055 — - a _ -z_
2
5.0E—6 “r -— - — _ __ _ _ ___________
L1 1"]
DOE-{>0 _. .l. u- -- I ........ I fl
0 500 1000 1500 2000
Pore width (um)

Figure 39. Comparison of motilities in pom with consumption and no consumption.

Small changes in the cell density profile in Figure 40 indicate that the interaction
from the pore walls doesn’t change very much in smaller pore sizes in the presence of a
large gradient. Three thousand cells responding to gradients in methylaspartate with a
chemotactic sensitivity of 4.1E—4 cmZ/s at 6 minutes were used for these comparisons.
Since the gradient is equally defined in the simulations, the enhanced migration through
pores must result from the combination of a larger gradient and a lower consumable

chemoattractant concentration.

63

 

 

— -o- - 1000 micron
—- 500 micron
—B— 300 micron
------>< 200 micron
—A—— 100 micron

x, = 4.154 (cmzls)

 
 
 
    
  

 

 

 

Relative cell density

 

 

 

Chamber position (cm)
Figure 40. Cell density profile for decaying step gradient in different pore widths.
6.2 MULTISLIT BARRIER
Accelerated migration was shown to happen through tubing connectors with a

pore width of 0.4 mm, both glass beads and sand also showed accelerated migration. To
explain the results of the experiment with cells moving faster through the porous material
then through the bulk phase, simulations were set up to model the experiment. A barrier
of parallel horizontal slit pores were set up on the right side of the simulation space with
the cells initiated at the center. These simulations will be used to explain the results of
the experiment in Figure 34. The porosity in the multislit simulations had a small effect
on the penetration ratio when the pore width was kept constant and only the solid region
width was increased. Simulations with larger solid regions azigure 41) did have lower

penetration ratios then the smaller solid region (Figure 42) simulations.

6000
4000
2000
Yunn) 0

-2000

-4000

 

-6000
-8000 4000 0 4000 8000
qun)

Figure 41. Migration into multiple loo-micron slit pores with a porosity of 33%.

6000
4000
2000
mm 0
-2000

-4000

 

-6000
-8000 -4000 0 4000 8000
qun)

Figure 42. Migration into multiple loo-micron slit pores with a porosity of 50%.
We have established that once a sufficient number of cells are in pores they can
collectively create large localized gradients causing accelerated migration in pores. Cells
outside of a porous barrier moving with only random motility will be greatly impeded by

the barrier and will barely penetrate into the pores. What causes the cells to be drawn

65

into the pores? Logic would state that there would have to be some kind of directed
migration towards the pore openings.

It was noted that a lot of bubbles formed in and around the porous media when the
cells migrated through them. Pseudomonas KC favorably consumes nitrate until it is
depleted then it metabolizes nitrite, which is a weaker chemoattractant. The metabolism
of nitrite results in the formation of nitrogen gas, which is most likely the cause of the gas
formation. The absence of bubbles in the rest of the bulk phase indicates that the
concentration of nitrate or nitrate is more plentiful then in the pores causing less
chemotaxis in the bulk phase. Low concentrations of multiple chemoattractants in series
may also cause accelerated migration in pores. Multiple gradients give rise to multiple
chemotactic waves of cells, which would increase the number of cells that could
chemotactically enhance their migration through the pores resulting in more cells on the
other side of the pores or further from their initial location. Our model currently can only
simulate one chemoattractant so enhanced gradients of a single chemoattractant was
studied.

The hypothesis that maybe there was additional nitrate in the sand core sample
that is drawing the cells in was tested by starting out with a higher concentration of
chemoattractant in the slit pores. This may explain the experimental results for the soil
but not the glass beads or tubing connector. This scenario seems plausible since in the
capillary experiments, a gradient of chemoattractant diffusing out of the capillary draws
the cells in to it.

As for the glass beads and tubing connector, what kind of concentration in the

pores could cause accelerated migration? We know (from the tumbling probability

66

equation and results) that chemotaxis is enhanced at low chemoattractant concentration.
If the-concentration of chemoattractant is lower in the pores this could cause accelerated
migration compared to the bulk even when going against the gradient at the mouth of the
pore.

What would cause the chemoattractant in the porous material to have a lower
concentration then the bulk concentration? Reviewing Laura Boom’s research notebook
indicated that the experiment that showed accelerated migration through the tubing
connector was wet when it was placed in the agar. The water that was in the connector
likely diluted the concentration of the chemoattractant in the connector or caused a dilute
chemoattractant film against the agar for the cells to swim through. Mark Widman’s
proposed explanation was the bacteria could swim faster due to reduced resistance
(1997).

These explanations seem reasonable since when the experiment was repeated with
a dry connector, accelerated migration was not noticed. The lab notebook does not
indicate if the beads and sand were dry after they were autoclaved. Residual water or
condensation may have caused the chemoattractant to be diluted in the porous media.

The accelerated migration could also be explained without having a diluted
chemoattractant in the porous media. The experiments showing accelerated migration in
the porous media lasted between 63 to 91 hours with a porous media thickness of 10.5
mm to 15.5 mm. It takes 2 hours for the number of bacteria in pores to double, which is
twice as long compared to growth in the liquid phase (Shanna, 1993). Even if a few cells

enter the porous media, the number of bacteria would double 45 times producing a lot of

67

bacteria in the pores. Results from the slit pore simulation indicate that a large amount of
cells in a pore will cause accelerated migration.

The results of our simulations indicate that increased swimming velocity or cell
growth and division is not required to accelerate migration through pores but may
enhance it. Modulation of the turn angle distribution caused by confinement of pore
walls is not needed either.

6.2.1 Non-reflective barrier

The first simulations were done with a uniform concentration of 0.1mM aspartate
and non-reflective boundary conditions when a cell encounters a wall. The consumption
rate used in these simulations is 0.1 5'1 cells". Different pore widths and slit wall
thicknesses were tried with different consumption rates to see if the cells would move
into the pores and move farther compared to the cells swimming in the opposite direction
in the bulk phase. In the simulations of cells entering multislit pores with widths of 900
microns in Figure 43 showed a penetration ratio of 2.06. The penetration ratios obtained
from experiments reported in Laura Boom’s lab notebook ranged from 1.2 to 1.6. The
radial cell density plot in Figure 44 shows a narrow chemotactic wave at its farthest point
and two other high cell concentration peaks were located at the mouth of the pore
entrance. Only one cell concentration peak developed for migration into the bulk area.

The concentration gradient is steeper in the porous area then in the bulk.

68

8000
6000
4000
2000
Yam) 0
-2000
.4000
-6000

-8000
-10000 -5000 0 5000 10000 15000
Xmm)

Figure 43. Cell migration into multislits with a uniform initial concentration.

 

 

 

 

 

 

 

 

 

 

9E+3 012
Cell density L—

8.E+3 -— * "~*~* — -CmmaMCMt——* P 0 1
1§7Ew—» — -*‘*~* ' A
g 6Ed3—r 'El§
., .5 e,
35.n3 E“
g 4E+3 - g g
==15E+3-~ g 5
<3 2Ea3- - 2% g

LE+3- 8

(1E+0

-18000 -12000 -6000 0 6000 12000 18000

Radial position (p, m)

Figure 44. Radial cell density profile and concentration profile with a uniform initial concentration.

69

Simulations with higher concentration of chemoattractant in the slit pores actually
show a decreased penetration into the pores. Although increased penetration occurs
when there is no chemoattractant initially in the bulk phase and it diffuses out of the
pores, which is the basis for the capillary experiments. These results may be explained
by the receptors being saturated in higher concentration in the bulk phase inhibiting
chemotaxis. The simulations with no chemoattractant initially in the pores resulted in a
lower motility into the pores indicating that a critical amount of chemoattractant is
needed in the pores to create an adequate gradient.

At the lower concentration in the pore, a gradient will have a larger effect on the
tumbling probability as equation 1 predicts leading to the increased motility through the
pore that we see in the simulations. Simulations in multiple parallel pores were done
with pore concentrations at half to a fifth of the bulk concentration. The solid pore walls
had a thickness of 100 microns. Pore openings ranging from 50 microns to 1800 microns
all showed an increase in penetration at 30 minutes compared to 15 minutes. A lower
chemoattractant concentration in the pores resulted in further penetration into the pores,
while a higher chemoattractant concentration decreased the penetration. Cells will swim
into the pore by random motility even though there is a gradient decreasing into the
mouth of the pore. Once they are inside they will create a gradient by consuming the
chemoattractant that will bias their migration through the pore. Figure 45 compares the
results of the simulations of cells entering multislit pores with widths of 900 microns
giving a penetration ratio of 2.24, which is 8 % larger then the in the simulations of
uniform initial chemoattractant concentration. The radial cell density plot in Figure 46

shows a smaller leading edge cell density peak and a large peak at the pore entrance. The

70

concentration profile was shallower due to diffusion of chemoattractant from the opposite

side of the pore openings.

Y(um) 0

8000
6000
4000
2000

-2000
-4000
-6000

-8000
400(1) -5000 0 5000 10000 15000

X<um)

 

Figure 45. Cell migration into multislits with a lower initial concentration in the pom.

 

 

   

 

 

 

1.6E+4 . 0.12
Celldensrty [
A 14E” “7"" i T F 'Tmifi - - -Chemoattractant 0_]
N: 1.2E+4 ee ,, e- , e V 1-, - , r“~~r~ A
31013444777777. 7| i, ,, A 77(770'08 EE
5 g c
i 8.0E+3e ,,,L,,,, k~W E g
33 '5
'3 6.0E+3— ,22 r-m e g
= F g :
<3 4.0E+3 - ._,L i-fiw . 3 g
D o
2.0E+3 ~— Lam, - o
0.0E+0 . T
-18000 -12000 -6000 0 6000 12000 18000

Radial position (p m)

Figure 46. Radial cell density profile and concentration profile with a lower initial concentration in

the pores.

71

The motility and penetration ratio both peak at a pore size of around 700 microns
and decreases at smaller pore widths and at larger pore widths as shown in Figure 47 and
Figure 48. The simulations with lower chemoattractant concentration in the pores always

resulted in a higher value.

 

 

Penetration ratio

 

 

 

 

 

0.5 ~ -~ , , —¥k-- — ,__,_L___
+0.1mMinpores
0 - -<> - 0.02mMinpores
0 500 1000 1500 2000

pore size (um)
Figure 47 . Effect of chemoattractant concentration on penetration ratio.

6.0E—5

 

5.0E-5 -.

4.0E-5 ~

 

(cmzls)

5. 3013-5 ..

2 2013-5 W

ouh

 

 

 

 

 

 

 

I'OE'S i' 8‘ i C i i“ “ " 7' +0.1mMinpores
- - 0.02mM'
0.01~:+0 . t '0 e "1me
0 500 1000 1500 2000
pore size (pm)

Figure 48. Effect of chemoattractant concentration on the motility coefficient.

72

Simulations were done using shorter pore lengths with a lOO-micron pore width to
observe the behavior when the cells exit a porous region. The acceleration of migration
continues when the cells are in the pores, but when the porous region stops the migration
is no longer enhanced. In Figure 49 the cells stay in the area against the outside of the
pores and are more likely to swim back into the pores when the concentration is lower in

the pores.

6000
4000
2000
YOU“) 0
-2000

-4000

 

-6000
-8000 -4000 0 4000 8000
qun)

Figure 49. Enhanced migration with initial lower concentration in multiple slit pores.

In Figure 50 with uniform initial chemoattractant concentration when the cells

exit the pores they fan out in a semicircular pattern and a chemotactic wave forms.

73

6000
4000
2000
YUM") 0
-2000

-4000

 

-6000
-8000 -4000 0 4000 8000
X(p.m)

Figure 50. Enhanced migration with initial uniform concentration in multiple slit pores.

In the short pores where the cells swim out the opposite end the highest

penetration ratio is when the concentration is the same in the pores as in the bulk

solution. The penetration ratio decreases at both higher and lower chemoattractant

concentration in the pores as shown in Figure 51, while the motility continues to increase

with decreasing concentration for the range of concentration simulated in Figure 52.

Penetration ratio

1.6
1.4 .-
1.2 .

1 -,
0.8 , , ~~~~~eve—envy , W, 4* ~~'~———~
0.6 we” ,, , we 777744 ,2 «447‘? ,,fl.,,i___

 

 

 

 

 

.0
h
1

 

0.2 77 7 7777 777777 77777v7
0 -9— Short 4,500 micron pores

0 0.05 0.1 0.15 0.2
Chennattractant concentration in pores (mM)

 

 

Figure 51. Comparison of the penetration ratio on short and long multislits.

74

3.0E—5

2.515-» H K

2.0E-5— 477 .7 -_ - LL

 

 

 

 

Motility (cmzls)
ii
1
l
l
l
i
l

 

 

 

 

 

 

LOB-5 7 —~ ~77 7~ A A- ,___. .,,_ __ __. -.._--__-._-._____

5-05'6 ‘ h a ”*7 7* -* 1 +Long18,000 nieronpores

0.0E+0 1 l-e—Short 4,500 rmcron pores .
0 0.05 0.1 0.15 0.2

Chemoattractant concentration in pores
Figure 52. Comparison of the motility coefficient on short and long multislits.
The effect of porosity of the porous region was determined by simulating 100-
micron pore widths with various sizes of solid pore wall between them. The penetration
ratio seemed to decrease then level off at higher porosity shown in Figure 53, while the

motility increased and then leveled off.

 

3013-5 3

2513.5 —- A _,__ . _ 2.5

1
N

2013-5 7 —7-— - —-—~———94 _ e- _ _M..._-L ___*-___,.2

 

Motility (cmzls)
I;
F.“
U!
l
1
1
j l
. I 1
? 3
O
>
.1
3.
Penetration ratio

1.0E-5 77 7777 7 ~77 — 4 —~ ~ m ~ ..____

1
fl

 

 

 

 

 

 

—o—O.l mM Motility
5.0E-6 7 7 7 7 77 e new”. - -O- - 0.02mMMotflity .. 0.5
—a—0.1 mM Penetration ratio
- -<>- - 0.02 mM Penetration ratio
0.0E+0 1 I 1 O
0 20 4O 60 80

Porosity (%)

Figure 53. The effect of porosity on the penetration ratio and motility in multislits.

75

These simulations only use one chemoattractant so when it is diminished biased
migration is inhibited. In experiments and natural conditions, bacteria are able to detect
multiple chemoattractant gradients and create new ones by metabolism of a primary
chemoattractant (i.e. nitrate to nitrite to N2).

6.2.2 Reflective barrier

For the case of reflective pore walls, the motility through the pores from our
simulation is higher than through the bulk when in the presence of a consumable
chemoattractant. In the cases where the chemoattractant is not consumable, the motility
through the pores corresponds with the motility in the bulk since movement in not denied
when it collides with a wall it is just repositioned. The reflective pores always resulted in
further penetration into the pores then the non-reflective simulations. The cell density
profiles for ZOO-micron pores in Figure 54 show that cells migrate farther into the pores

when the pore walls are reflective.

 

l .2E+4
— Reflective wall

1.0E+4 77777 777— _ ~ — —~ - ~~7Am_————_—___. - - - man-reflective wall

 

 

 

8.0E+3 e7——~ew-—~

 

 

6.0E-1-3 «4 “”24“”

 

 

4.0E-1-3 - ,,_-._~_.__________

Cell density (Cells/cmz)

 

2.0E+3 4* A. — - -- - ,

 

 

 

 

 

0.0E+0 r i
~18000 -12000 -6000 0 6000 12000 18000

Radial position (pm)

Figure 54. Comparison of reflective and non-reflective walls of a Zoo-micron multislit.

76

The mode found in natural environments is somewhere between the reflective and
non-reflective case. The type of soil and the adsorption properties of the bacteria would
determine the amount of retention of the bacteria in the pores. Figure 55 shows the cell
snap shot of cells moving through multislit pores with 900-micron width with reflective
pore walls. The penetration ratio is 2.33 which compares to the non—reflective

simulations where the penetration ratio was 2.06, yielding a difference of 11.6 %
6000

4000

2000

Y(llm) 0
-2000

-4000

 

-6000
-6000 -1000 4000 9000 14000

Xmm)

Figure 55. Simulation results of reflective 900-micron multislit pores.
The motility into the pores compared to bulk motility increased with time for pore
sizes from 100 micron and larger for all the concentration ratios simulated with 0.02 mM,
0.1 mM, and 0.2 mM in the pores and a uniform concentration of] mM every where else.
The simulation using 50—micron pore size and concentration in the pores a fifth (0.02
mM) of the bulk initially showed increased penetration ratio at 15 minutes, but at 30
minutes the penetration ratio decreased but was still larger then in the bulk. Since the

boundaries are reflective the decrease is not more resistance with the pore wall. This

77

enhancement of the motility is only temporary and only lasts as long as the cells are in
the pore. Once the cells reach the end of the pore the concentration increases and the
biased migration decreases. In the absence of growth, there are relatively few cells at the
end of the pore, which takes longer to produce sufficient local gradients and a lower
concentration.

The lOO-micron pore with equal Figure 56 had the largest penetration ratio of
3.29. For the reflective wall boundaries, the case with higher and lower concentration in
the pores resulted in a lower penetration ratio. The simulation with 0.02 mM
chemoattractant concentration in the pores is a lot lower then the bulk and results in a
higher penetration ratio at 15 minutes and a lower penetration ratio at 30 minutes

compared to the constant chemoattractant concentration.

6000 -
4000 - «’
2000 - '
Y(um) 0 n

-2000 e

 

-4000 - I? .

 

-6000 —
-6000 -2000 2000 6000 10000 14000 18000
X(um)

Figure 56. Simulation result of loo-micron pore width multislit with reflective walls.

Figure 57 and Figure 58 show the relation between the 2D and 3D motion on the
simulations in the slit pores. The 3D simulations showed that the motility is lower then

the 2D simulations for different pore sizes for both the reflective and non-reflective

78

 

migration through the pores with a difference of about 607 . The penetration ratio for the

2D and 3D simulations correlated pretty well for both the reflective and non—reflective

simulations with a difference of around 6%.

 

3.5

'1 :d E 2323.36.16. Massage: 2:

~\\\~\\\\\\ss\\\s«\s\‘\\-
\~\\\\\\\\\\\\\\\\\\\\\\\\
\\\\\ \ ~\\ «s u
\\\ \ \
Rt

~\\\\\s\\\\\\~\
\ \\\\\

..si. (\ss....ss.,.
s .s.................................. .e.
‘ \\:‘\

. \\
§s\ \\

\\\\\\\s\\\fl\\\\\\\ \
\ \ \ \ ‘\
=\\:‘\\\\s\\:\:=:::\

(tsSss

\\ \
SSRS?
“$33:

=l.:.:

g$§

\\ \ \\ \ \\ \\\\\\\\\\\\\s \\\\\
\\=\=~\§=\=ss=\\\\\\\\\\\\s\:\\\\\\\

\ \\ \ \\\\s \\\\\\\ \\§\
‘§‘§\\2§eteseeeks\\ess\¢§\s

s
\\~\\\\\\\\\\\
\\\\\\\\\\ \
::::::::2 ~

5 ‘ \x‘ 1
\x\\ “\\ \\\\\\‘

s
s
\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
\\\\\s\\\\\\\\\\~\
\\\ \ \\
\\
s

s

i\\\\£\'\\\'\\'\\'\' \\'
..
‘ttt

\\\s~ \\
\‘a\\\\\

b

 

 

 

0.1 leI 0.02 mM 0.1 mM 0.02 mM
200 mbron 200 micron 900 micron 900 micron

Chemoattractant concentration in pores (mM)

Figure 57. Comparison of 2D and 3D motion on penetration ratio.

 

7.0E-5
fl I 3d 3 2d V/é Reflective 3d E Reflective 2

6.0E-5 e W

vrnm
nun

33:36::
.........
\ ititzztztiii:tz
.....,...........

HER?
\\\\\

Motility (cmzls)

 

 

 

0.1 mM 0.02 mM 0.1 mM 0.02 mM

200 micron 200 micron 900 micron 900 micron

Chemoattractant concentration in pores (mM)
Figure 58. Comparison of 2D and 3D motion on motility coefficient.

79

6.3 HETEROGENEOUS
Three different effective shapes of particles were simulated. The first is circle
disks (Figure 59) with 2D motion, which is completely packed and closes off migration at

22% porosity.

 

Figure 59. Cell trajectory through packed bed of 2D disks

The other two shapes use 3D motion and are a cylinder and a sphere. The highest
porosity possible using spheres on our square lattice grid is 48%, which corresponds to a
face centered configuration. In the 3D motion simulations the z-axis uses a periodic
boundary so migration is through a bed of long (Figure 60) or strings of spherical beads
(Figure 61) respectively. When the cylinders are completely packed 22% porosity, the

cell can still move in the z direction but is hindered in the x and y direction.

80

 

Figure 60. Cell trajectory through packed bed of cylinders

 

Figure 61. Cell trajectory through a packed bed of spheres.

81

Where as the spherical string of beads can be completely packed and the cells can
move through the interstitial spaces as in a packed bed of spheres. Figure 62 compares
the motility coefficient of 10000 cells starting out in the center of a bed of the three
different shape particles with SOC-micron diameters. The consumption rate used in all the
heterogeneous porous simulations is 0.1 s’I cells". The migration is faster through the 2D
disks using 2D motion then through the 3D spheres and cylinders using 3D motion,
which is expected in the presence of a 2D chemoattractant gradient. The motility for all

three of the particle shapes decreased dramatically at about the same porosity of 50 %.

 

 

1.6E-5 . ~-—— * —

 

 

  

 

Motility (cmzls)
m t—t
o N
in in
0\ UI

l
l
l
l
i
I

 

 

 

 

 

 

4013-6 4.- _ _ 7 e -. _ -I- 3D Cylinder 500
/ J -o- 2D Disk 500
0.0E+O , 1 I -&- 3D §pmm 500
0 20 40 6O 80 100
Porosity (%)

Figure 62. Motility of different shape 500-micron particles
6.3.1 Migration into a strip of porous media
Randomly placed particles were added to the simulation grid to simulate
migration into porous media. Simulations using a strip of porous media on the right side
of the simulation space were used to correlate with the results of the multislit simulations

as well as previous experimental results. Figure 63 shows results of a simulation of

82

10000 cells starting at the center and moving into a strip of 400-micron particles with 3D
motion. The black dots are the cells and the gray circles are the particles with a porosity

of 57%.

4000
3000
2000
1000
YOU“) 0
-1000
-2000

—3 000

 

-4 000
4000 -2000 0 2000 4000
X(tlm)
Figure 63. Simulation results showing cells migrating into porous region.

The porosity was calculated from the area not taken up by the particles. The
average porosity of the Schoolcraft bioremediation site sand is 40.2 % (Criddle, 1997).
Experimental tortuosities for particle sizes of 100 to 800 ranged from 4 to .5 (Barton,
1995). The tortuosity increases as the particle size decreases because of the smaller
effective pore diameter and the curvature of the pore is larger and causes a collision with
the pore wall before its basal run length. When cells react to a chemoattractant gradient,
the run length is increased effectively increasing the tortuosity (Barton, 1995). We

observed that the motility through the pores decreased with decreasing particle size.

83

Experimental results indicate that the particle size has an effect on the motility
independent of the porosity (Duffy, 1995).

In Figure 64 the cells were started in a line along the edge of the porous region
instead of a point to see how the initial location of the cells affected migration into porous
regions. The cell movement and chemoattractant diffusion is periodic in these
simulations. Two bands of high cell density form, one migrates into the porous region

and the other migrates away from the porous region.

10000 7

5000 7

 

Yum) 0 7

-5000 -

 

 

—10000 7
-10000 -5000 0 5000 10000

xmm)
Figure 64. Simulation result of cells starting out in a line along the porous region.

Starting the cells out in a line along the edge of the particles versus starting them
out at a point resulted in a lower motility and penetration ratio as seen from the cell
density profile in Figure 65. The rest of the simulation results are done with all the cells

staring out at a point.

84

4'OE+31— ~7-7—7 77 ._ 7

 

— Initial line of cells
3,554.3 7, A_L__h - 3A -G—InitialJ)oint of cells —,
h

 

 

 

3.0E+3 7
2.5E+3 7 777777-
201348 77 7777
1.5E+3 77777
1.0E+3 7 77
5.0E+2 j - ‘

Cell density (Cells/cmz)

 

 

 

 

-10(X)0

Radial distance (um)
Figure 65. Cell density profile for initial cell location of point versus line.

Simulations were done using particle diameters ranging from 200 microns to 1600
microns. The motility coefficient for the cells moving into the pores and the cells moving
away from the particles was calculated separately. Using a porosity of 65% at 2 hours,
the motility coefficient changes very little over the range of particle sizes from 200 to

1600 microns as shown in Figure 66.

 

 

 

 

 

 

 

1.6E—5

1.4E-5

7; 1213-5 4““ —'
E m... 1183353315123:
we. manager

 

 

 

 

§ 6.036 — W

2 4.036 .- M 5m.
2.036 43 .m A e - e - - e - w..-“ __ __s
0.03+0 . i .

0 500 1000 1500 2000
Particle diameter (11m)

 

 

 

 

 

Figure 66. Motility results for different particle diameter in strip.

85

While at a porosity of 40 %, the motility coefficient for the cells moving into the
particles decrease at smaller diameter particles and the motility coefficient of the cells
away from the particles increases slightly since more cells move away from the particles
due to obstruction in the opposite direction by the particles. In Figure 67 the penetration
ratio seemed to peak before decreasing more at smaller diameter particles. At a porosity
of 65% the peak in penetration ratio was around a particle size of 500 microns, while a
porosity of 40 % peaked around 1000 microns.

1.00

0.90 ---.____
0.80 .1 _ __
0.70 at. -- -

 

 

 

 

 

 

 

 

 

 

 

 

e
a: _
g .1..- 127 __ ._~_ 3
.3 0.50 F b ' F #’—‘_flh“wwmfign__fl
g 0.40 777 7 7 777777777777777w7—7—
E 0.30 777 7 — 77——7—
0'20 ““5 5’ 5 “WW—”‘7 +65% Porostiy penetration
0‘10 2‘ —i v - ‘"_ —__—#— -B—40% Porosity penetration
0.00 1 1 1 1 1 1 1
0 200 400 600 800 1000 1200 1400 1600
Particle diamter (pm)

Figure 67 . Penetration ratio result for different particle diameter in strip.

The penetration ratio and the motility coefficient of cells swimming into the
particles both decreased at lower porosity, while the motility coefficient for the cells

swimming away from the particles increased at lower porosity as seen in Figure 68.

86

 

1.20 1.4E-5

 

 

 

 

 

 

 

 

 

 

1.00 g g L, f k __g____, 1.235
.2 , A
a 0.80 _ __ ,_ ik,i vs a, _‘ .~ IDES N£
E E
1: 77 8.0E-6 3
'3 0.60 7 777 7 7