|UFDC Home||myUFDC Home | Help|
This item has the following downloads:
SEROTYPE REPLACEMENT OF VERTICALLY TRANSMITTED DISEASES
THROUGH PERFECT VACCINATION
DOUGLAS H. THOMASEY
A THESIS PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
MASTER OF SCIENCE
UNIVERSITY OF FLORIDA
2007 Douglas H. Thomasey
I dedicate this graduate thesis to my parents Rich and Laura, my brothers Bryan and
Kevin, and my girlfriend Ashley. I could not have done this without your support.
I would like to thank my advisor, Dr. Maia Martcheva, and my committee members,
Dr. Sergei Pilyugin and Dr. Be( i v,:iii Bolker. This would not have been possible without
your insight and guidance.
TABLE OF CONTENTS
ACKNOW LEDGMENTS ................................. 4
LIST OF TABLES ....................... ............. 7
LIST OF FIGURES .................................... 8
ABSTRACT . . . . . . . . . . 9
1 INTRODUCTION ................................ 10
1.1 Vertically and Horizontally Transmitted Diseases ............. 10
1.1.1 Eradication of Horizontally Transmitted Disease through Vaccination 11
1.1.2 Vertically Transmitted Disease with Vaccination . .... 12
1.1.3 Importance of the Reproduction Numbers . . 13
1.2 Multi-strain Diseases and Vaccination .................. .. 13
1.2.1 Mechanisms of Strain Replacement . . . 14
1.2.2 Similarities between the Known Mechanisms of Strain Replacement 15
1.3 Hypothesis ................ .............. .. 16
2 VERTICAL TRANSMISSION MODEL WITHOUT VACCINATION ...... 17
2.1 The General Model . . .... ... . . 18
2.2 Model 1: A Model without Healthy Births from Infected Individuals . 20
2.2.1 Disease-Free Equilibrium ............... .. .. 20
2.2.2 Reproduction Numbers .............. ... ... .. .. 21
2.2.3 Strain 1 and Strain 2 Equilibrium ............... .. .. 24
2.2.4 Invasion Numbers ......... . . . .... 28
2.2.5 Expressing the Invasion Reproduction Numbers as Functions of the
Reproduction Numbers ................ ... ... .. .. 31
2.2.6 Parametric Plot . . . .. ... ......... 32
2.3 Model 2: A Model with Healthy Births from Infected Individuals . 40
2.3.1 Equilibrium and Reproduction Numbers . . 40
2.3.2 Strain 1 and Strain 2 Equilibrium ............... .. 41
2.3.3 Invasion Numbers ............ . . ... 43
2.3.4 Simulation Plots ............... ......... .. 44
3 VERTICAL TRANSMISSION MODEL WITH VACCINATION . ... 46
3.1 Model with Vaccination . . . .... . .. 46
3.2 Model 3: A Vaccination Model without Healthy Births from Both Infected
Individuals ................ ............. .. 47
3.3 Model 4: A Vaccination Model with Healthy Births from Both Infected
Individuals ................. ............ .. .. 52
4 DISCU SSION . . . . . . . .. 54
REFERENCES .............. ........... ... ... 57
BIOGRAPHICAL SKETCH .................... . .59
LIST OF TABLES
2-1 Parameter meanings for the model ............. ... ..... 17
4-1 Trade-off mechanism overview ............... ........ .. 54
LIST OF FIGURES
2-1 Flow chart of the model without vaccination . . .....
2-2 Parametric plot of the invasion numbers in Model 1 . . .
2-3 Shaded parametric plot of Model 1..... . . .....
2-4 Parametric plot of Model 1 highlighting a specific area of suboptimal reproduction
number dominance ....................... . . .....
2-5 Simulation plot showing a strain with lower reproduction number being dominant
2-6 Parametric plot of Model 1 when we do not get the situation of the strain with
suboptimal reproduction number being dominant . . .
2-7 Existence of a unique solution using the Intermediate Value Theorem . .
2-8 No solution for J* using the the Intermediate Value Theorem . . .
2-9 Simulation of Model 2 when the strain with larger reproduction number
3-1 Flow chart
of Model 2 when the strain with suboptimal reproduction number
of the model with vaccination. . . .
and simulation plots of Model 3 when ) = 0 . . .
and simulation plots of Model 3 when = 0.15 . . .
and simulation plots of Model 3 when y = 0.28 . . .
of Model 4 when ) = 0......... . . .
of Model 4 when i = 0.2........ . . .
of M odel 4 when p = 0.23 .............. . . .
Abstract of Thesis Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Master of Science
SEROTYPE REPLACEMENT OF VERTICALLY TRANSMITTED DISEASES
THROUGH PERFECT VACCINATION
Douglas H. Thomasey
('C! i: Maia Marcheva
There have been many confirmed mechanisms which lead to strain serotypee)
replacement. This replacement is the result of the dominant strain within a population
being replaced by another strain, that initially was the non-dominant strain. Differential
effectiveness of the vaccine is the widely accepted and the most important mechanism
which leads to this replacement effect. However, in this study, we will inspect a concept
known as "perfect" vaccination. Perfect vaccines are vaccines that are assumed to be
one-hundred percent effective against both individual strains equally. It has already
been proved that "perfect" vaccination, along with a trade-off mechanism, such as
super-infection or co-infection, lead to strain replacement. But now, the focus will be to
examine if this same effect is possible in regards to a model where vertical transmission is
the trade-off mechanism.
Until recent times, biomathematics was not considered to be one of the popular fields
in mathematics. However, it is becoming more and more renowned to mathematicians
across the world. There are many reasons why this type of mathematics is important, and
within this research we will demonstrate why. Mathematical modeling gives us the ability
to foresee how a certain problem may evolve over time. Different detection, prevention,
and control strategies can be tested and evaluated before they are implemented . In
respect to diseases, this method proves to be a very powerful tool in the aid of preventing
a disease spread.
1.1 Vertically and Horizontally Transmitted Diseases
Infectious diseases are categorized in a number of different v--~i-. These categories
may vary from which part of the body the disease affects, its symptoms, what causes
the disease, to whether or not the disease is curable. However, another important
classification of diseases comes through their means of transmission. There are two
V--IV- of transmission, horizontal and vertical. Horizontal transmission is defined to be the
passage of the infection from one host individual to another . This can occur by direct
physical contact and by the inhalation or ingestion of infective material .
The category of diseases that interests us and that will be examined throughout this
study, are the vertically transmitted diseases. These are the types of diseases that use a
direct transfer of a disease from an infective parent to an unborn or newly born offspring
. This can occur transplacental (or through infected seeds when dealing with plants)
and by contact with blood or other bodily fluids, resulting from sexual activity and needle
sharing , . Horizontal and vertical transmissions differ greatly in the way they are
spread and how they affect humans. Compared with the pure horizontal transmission case,
we do not have many results for vertically transmitted diseases in structured populations
. Examples of vertically transmitted diseases include chronic illnesses such as acquired
immunodeficiency syndrome (AIDS), human T-cell leukemia virus (HTLV), Hepatitis B,
C'! ;, I-' disease, rubella virus, and herpes simplex virus , . Much of what we found
may be very helpful in understanding how these diseases function, as well as what may be
the best methods to go about preventing their spread.
1.1.1 Eradication of Horizontally Transmitted Disease through Vaccination
Vaccination p1 i' a crucial role in controlling the spread of disease. Without them,
diseases are able to spread uncontrollably and have the potential to become a pandemic.
An example of how vaccines have helped in the past can be seen through the eradication
of smallpox. Even though this is a horizontally transmitted disease, it allows us to
understand easily how vaccination can bring a disease to eradication.
In endemic form and in waves of epidemics, smallpox killed and disfigured millions
of people throughout the world . In epidemic years alone, between 1868 and 1885,
smallpox killed an estimated ten percent of the population of British India . However,
due to the unselfish action of the The World Health Organization (who initiated a
program in 1967 for the global eradication of smallpox and achieved it in October 1977),
we saw the last person to acquire naturally occurring smallpox in the world recover from
this disease in Somalia, Africa . Even though this was a great accomplishment for
the human race, some people understood the power of smallpox and considered using
developed strains as a weapon for mass destruction .
Recently, the United States recognized the threat of the avian flu, whose more
common name is the "bird flu". The main reason for this fear has been brought about
by the confirmed outbreaks with human fatalities in nine different Asian countries since
2003 . In addition to infecting humans, there have also been numerous outbreaks in
poultry farms. Some of these reports came from Hong Kong in 1997, as well as South
Korea in 2003 , . With no vaccine for the avian flu, if it becomes human-to-human
transmittable and spreads among humans, it will be hard to estimate how many people
would die. For this reason, in 2006 the people in the United States recognized this
problem was real. For now, there may be more emphasis towards problems we could face
with our food source, especially poultry farms, than humans. But nonetheless, the idea
that humans could get a disease with no vaccination available is now evident to everyone.
As we look to the future, vaccines will prove to be more and more valuable in the
fight against disease. The pathogens that cause certain diseases, such as the flu virus,
seem to find a way to mutate and develop into new strains. For this reason, the idea of
evolution as a whole can be considered a natural strain improvement program . With
this in mind, the more we know about developing vaccines, the better off we will be in
producing new vaccines. Scientists have determined an effective way of using gene cloning
to produce novel vaccines . When dealing with the gene cloning process, we do not
include specific agents in a vaccine. Instead, we have found a way to include more versatile
genes into a vaccine in order to vaccinate against a particular disease more effectively.
This idea can then be used to include genes from several different infectious agents within
one vaccine, thus raising the possibility of using a single vaccination to confer protection
against a number of diseases simultaneously . Without knowing, we may be vaccinating
against diseases in which we did not know were present.
1.1.2 Vertically Transmitted Disease with Vaccination
The obvious reason for studying vertically transmitted diseases is to find a way
to cure them. From the previous section, we saw that vaccination offered the ability
to eradicate smallpox. But, this raises the question if the same can be done with the
vertically transmitted diseases? Hepatitis B and rubella are examples of vertically
transmitted diseases which are vaccine preventable . On the other hand, diseases
such as AIDS, chagas, and gonorrhea are vertically transmitted, but there are no vaccines
for them yet . Even if vaccination was capable of eradicating these diseases, it is left
to humans to develop a vaccine that would be effective. In addition to this, even with
a vaccine against these diseases, we may not be certain this would solve the problem.
Biomathematics and mathematical modeling may allow us to make predictions before
spending a lot of time and money for something that will fail.
1.1.3 Importance of the Reproduction Numbers
The reason why vaccines are so powerful against eliminating a disease comes from
the attempt to reduce a particular diseases basic reproduction number below the critical
level of one . The reproduction number indicates the number of secondary infections
one individual will produce in its lifetime in an entirely susceptible population. With this
in mind, if it is possible for a vaccine to reduce the reproduction number below one, then
in most cases the disease will die off. The difficulty of developing these vaccines comes
from the fact that our immune system has a diversity of 108, which makes it very hard
to find the correct makeup of the vaccine . Hopefully, with the understanding of how
these disease react in certain situations, we will be able to find vaccines for them more
1.2 Multi-strain Diseases and Vaccination
When a particular disease only has one form, or strain, vaccination is a very powerful
tool to eradicate the disease. However, vaccination becomes much less efficient when a
disease has more than one strain, particularly if it has hundreds of different strains. Some
examples of diseases with more than one strain include Hepatitis B, Tuberculosis, as
well as forms of Diarrhea and food poisoning. On a higher scale, a disease that possesses
one of the largest number of different strains is that of influenza, which has hundreds of
different strains. Most notably are the strains causing the Spanish, Asian, and Hong Kong
Influenzas, as these strains, when first introduced to the human population, developed into
pandemics killing more than 20 million people .
The more strains a certain disease possesses, the harder it may be to vaccinate against
this disease. Not only do you have to vaccinate against the correct disease, but also the
right strain. If many strains circulate in a population and we vaccinate against the most
dominant strainss, the number of cases with this strain declines. In the absence of this
strain, the other strain, which was not dominant before, now may become the dominant
strain and bring with it a new set of problems or symptoms. This process is known as
population-level strain serotypee) replacement.
All types of serotype replacement require some type of mechanism for the replacement
to actually occur. The main mechanism of serotype replacement is through differential
effectiveness of the vaccine. This is the mechanism responsible for strain replacement
based on the fact that vaccination against a certain strain eliminates it; however, other
strains previously not as prevalent, and against which the vaccine is only partially efficient
or not efficient at all, now emerge and become the dominant . This leads us into the
main focus of this study and the fact that mathematical models -i-i-.- -1 that serotype
replacement may also occur through "perfect" vaccination.
Perfect vaccination is the process of vaccinating against all different strains of
a virus and the vaccinated individuals are one-hundred percent protected against all
strains . With perfect vaccination, serotype replacement may occur if some trade-off
(coexistence) mechanism is in place. A trade-off mechanism is a mechanism that allows
for the coexistence of two strains (within a host or on population level) . The reason
for this comes from the idea that different amounts of a perfect vaccine applied to a
population have different affects. No vaccination or even a small amount of vaccine
would not influence the dominant strain enough and it would still dominate. But as we
increase the levels, vaccination should be able to permit the two strains to coexist, hence
a trade-off mechanism would be necessary. However, strain replacement would not occur
until the vaccination level was raised enough to have the non-dominant strain take over
and become dominant. For this reason, we may assume that a trade-off mechanism is
necessary for perfect vaccination to be a mechanism for serotype replacement.
1.2.1 Mechanisms of Strain Replacement
Until now, it has been established that there are two trade-off mechanisms that can
facilitate strain replacement with perfect vaccination. Super-infection has proved to be
one such trade-off mechanism . Super-infection is the process where an individual
infected with one strain of a virus comes into contact with an individual infected with
the second strain, becomes infected with the second strain, which takes over the host
immediately . It has also been proved that perfect vaccination and the trade-off
mechanism co-infection also leads to serotype replacement . Co-infection is the process
where an individual infected with one strain of the virus comes into contact with an
individual infected with a second strain, not necessarily of the same virus, and becomes
infected with both strains at the same time .
1.2.2 Similarities between the Known Mechanisms of Strain Replacement
When looking at the two trade-off mechanisms that lead to strain replacement, we
notice there are some similarities between them. One of these similarities is that both
mechanisms have the property to allow a strain with a lower reproduction number to
exclude a strain with a higher reproduction number (in the absence of vaccination).
Even though in the absence of these mechanisms the strain with the higher reproduction
number will eliminate the strain with the lower reproduction number, these mechanisms
provide the means necessary to enable dominance of the strain with suboptimal
reproduction number. The reproduction numbers prove to be a very useful tool in
analyzing a particular system. Different aspects of a system including equilibria, stability,
as well as dominance, all can be examined through the use of the reproduction numbers
. In addition to these uses, the reproduction numbers also serve as the basis for finding
another meaningful expression known as the invasion numbers. An invasion reproduction
number gives the number of secondary infections that an individual infected with a strain
will produce in a population where the other strain is at equilibrium .
The second similarity between these two trade-off mechanisms is that they both are
the result of with-in host competitive interaction of the strains, or the interaction of one
strain of a virus in a particular host with another strain in the same host . With this
in mind, we can understand the importance of results acquired as they deal with a single
host. If it is known that the pathogen of a certain disease has more than one strain, all
which are vaccinated against, we can now understand the outcome this vaccination will
have on both strains within the single vaccinated host.
After understanding what has been done in the past through the introduction, we can
now present the main hypothesis which will be addressed in this study. The main question
to be addressed is what trade-off mechanisms lead to strain replacement with perfect
vaccination? The hypothesis which stems from this question is based on the importance
of the reproduction numbers. We expect that only trade-off mechanisms which lead to the
dominance of a strain with suboptimal reproduction number lead to (population level)
strain replacement with perfect vaccination.
To address this hypothesis, in the following chapters we will examine a few different
mathematical models and interpret them using different means, which include analytical
interpretations and simulations. The structure of the models come from a paper written
by Marc Lipsitch . The baseline model is a model with two strains and vertical
transmission. Vertical transmission pl .il- the role of a trade-off mechanism, and allows
for coexistence. We will examine the baseline model in C'! lpter 2. We will use two main
forms of a vertical transmission model, that is, we will break this model down into two
cases. First, where we do not have healthy births from infected individuals, and second
where I. !lrl !: births from the infected classes are permitted into the susceptible class.
These two models will be studied in C'! lpter 3, however vaccination will be included
within them. We utilize models without and with vaccination to get an idea of what may
occur when there is no vaccination class present, in contrast to the results we may find
when we actually have a vaccinated class. Finally, in Chapter 4, we will discuss all of the
findings and what they actually mean in connection to vertically transmitted diseases.
VERTICAL TRANSMISSION MODEL WITHOUT VACCINATION
In this chapter, we introduce a model with vertical and horizontal transmission
without vaccination. We consider both the general model and a plausible simplification.
However, through the different assumptions, the two versions of the model are mathematically
different. With this in mind, we will consider each case separately and analyze the
different effects each model produces. We will begin by introducing the general model.
We consider a population, with respect to time, of total size N(t), which is then
divided into three different classes. The number of susceptibles, uninfected with either of
the two strains, is denoted by S(t). The two infected classes are represented by I(t), for
those infected with strain one, and J(t), for those infected with strain two.
Table 2-1. Parameter meanings for the model
bx per capital birth rate into the susceptible class
bi per capital birth rate of infected newborns into infected class I
b2 per capital birth rate of infected newborns into infected class J
rl per capital birth rate of uninfected newborns from infected class I
r92 per capital birth rate uninfected newborns from infected class J
93 per capital birth rate of susceptible newborns from vaccinated class
(to be used in the vaccination model)
31 transmission rate of strain 1
32 transmission rate of strain 2
p per capital death rate
a disease-induced per capital death rate of infected class I
a2 disease-induced per capital death rate of infected class J
b per capital vaccination rate (to be used in the vaccination model)
In addition to the different classes, there are many different parameters included
in the model. (Table 2-1) is a reference for what each of the parameters within the
model corresponds to. From (Table 2-1), we can understand what influence most of the
parameters have on the model, such as the birth and death rates. However, it may be
helpful to explain a few of the other parameters in more detail. For instance, the rate
susceptible individuals become infected with strain one is given by the transmission rate,
31, of strain one. Symmetrically, 32 signifies the rate at which healthy individuals become
infected with strain two. The other two parameters that may not be so easily understood
are ac and a2. Typically, an individual infected with a particular virus has a higher death
rate than an individual in the susceptible class. Since we have assigned the natural death
rate to be p, al and a2 give us the ability to increase the death rates within the two
2.1 The General Model
With the notations above, the general model is as follows.
(bS + I + I/2J)(1 S I J) pS P1IS 02JS
I[bl(1 S I J) (p + a,) + S]
J[b2( S I J) ( + a2) + 2S]
The term (1 S
total population N(t)
I J) represents the carrying capacity of the system. Also, the
S(t) + I(t) + J(t) satisfies the following differential equation.
d = (bS + rl + T12J+ bil + b2J)(1 S- I J) ipN all a2J.
To understand how the general model works, we can represent it as a flow chart.
Figure 2-1. Flow chart of the model without vaccination.
Assumptions. In order for certain conclusions to be made, we must offer a few
assumptions in regards to the parameter values. From a biological standpoint, it is
reasonable to assume that the birth rate of new individuals from either of the infected
classes is less than that of the healthy population. Thus, we are lead to the following two
r91 + b < b, (2-3)
2 + b2 < b, (2-4)
Having these two assumptions then leads us to another important assumption through
using the total population equation. We can express (2-2) in terms of N as the following.
Sdt (q +b) I+( +b2)J](1 S -I -J) pN all J
< b1N(1 N) pN (Logistic Equation)
< bN iN
= (b. p)N (\! I1 Equation)
~ N() < Noe(bx-p)t (2-5)
From the preceding calculations, we now have the inequality (2-5), which gives us a
bound of the total population size at time t. It can be observed, from the two parameters
remaining in the equation, that we have two cases. The first case is where b, < p. With
this assumption we would have a decreasing exponential function, and N(t) -- 0 as
t -- oo. So the entire population would die out. This case is not very interesting because
we expect the population to persist in order to investigate the dynamics of the disease in
it. For this reason, we will make another assumption for the initial parameter values.
b, > p (2-6)
First, we will analyze a model where healthy newborns are not allowed from either
of the two infected classes. Then we will explore the same model, however now with the
assumption that healthy newborns from the two infected classes will be incorporated into
The two cases can be expressed in terms of initial parameter conditions by the
1.) There are no healthy newborns from either of the two infected classes, r11 = T12 = 0.
2.) There are healthy newborns from either or both infected classes, T11 / 0 and\or
/2 / 0.
2.2 Model 1: A Model without Healthy Births from Infected Individuals
With the first model of interest we assume r11 and 172 to be zero, which should be
noted before we start analyzing. For this model, we will go through the analysis in detail
to give the reader an understanding of the steps necessary to carry out full analysis. In the
later models we may just list certain expressions, assuming one can find the expressions in
the same manner.
2.2.1 Disease-Free Equilibrium
The first step in analyzing a model is to compute the the disease-free equilibrium
(DFE). This will give us a particular ordered triple, of the form (S*, 0, 0), to signify
the equilibrium level of the susceptible class when disease is not present. Observe
that "S*" represents a time-independent variable. It is a stationary value and also a
time-independent solution. Nonetheless, finding the DFE is normally the first step in the
analytical process as it offers an important substitution for the reproduction numbers. The
process of finding the DFE starts by setting the derivatives equal to zero in the original
system. For the DFE, we also have I* = 0 and = 0.
0 S*b, S*2 S*1
Clearly, the above system (2-7)
solve for S*.
0 = b2S*
*2 = bxS*
is trivially satisfied. Using the first equation we can
From this we obtain the DFE ordered triple, So (= ,, 0, ).
2.2.2 Reproduction Numbers
The next step in the computational process includes finding the jacobian matrix. This
is done by taking the derivative of each equation with respect to one variable per column.
For example, the first column of your matrix will have entries from the three equations
differentiated with respect to S. The second column will have the same orientation,
however the equations will be differentiated with respect to I. This process would then be
repeated for the number of equations that the model possessed. From system (2-2), we
have the following matrix.
( b(1-2S- I -J) -31-32J
-b2J + 32J
bi(1 S- 21- J) (p + ai) +3S
b2(1-S-I-2J)- (p+a2)+ 32S
By way of the jacobian matrix, we now want to compute the reproduction numbers
R1 and 72 of the system. In order to do so, we have to compute the jacobian at the DFE.
Here, it is necessary to evaluate the jacobian using the DFE ordered triple (S*, 0, 0) and
then find the eigenvalues of the resulting matrix.
b,(1 S*) bS* -bS* bS* bS* 02S*
0 bi(1 *) (p +a) + lS* 0
0 0 bo2( S*) ( a + 02 S*+
In recognizing that matrix (2-9) is upper triangular, the eigenvalues are the values of
the diagonal elements. We can label these values as A1, A2, and A3 respectively. But more
importantly, these eigenvalues give us the expressions needed to compute the reproduction
Mathematically, the reproduction numbers, when less than one, provide conditions
so that A1 < 0, A2 < 0, and A3 < 0. Thus, to compute the reproduction number for
strain one, or R1, we will use the eigenvalue from the column of matrix (2-9) which was
differentiated with respect to I. In this case, we will use A2, or entry (2, 2), from matrix
(2-9). The condition A2 < 0 is equivalent to the following condition.
bi(1 S*) (p+ a)+ + S* < 0
bi(l S*) + 1S* < (P + a1)
bi(1 S*) + 1S*
bi(1 S*) + -- S*
p + a,
We can compute the reproduction of strain two, 72, in the same manner. However,
in addition to finding the reproduction numbers, we can simplify them by substituting
parameter expressions for the terms involving the different class identifiers, such as S* and
(1 S*). This substitution comes from the system for the DFE (2-7). In having already
solved and simplified this expression (through finding the DFE), when we substitute the
expression for S* and (1 S*), we obtain the reproduction numbers.
bl + (b p) b2P +32(b, (p)
R = 2 = (2-10)
b.(p + a,) b.(p + a2)
After finding the two reproduction numbers, we notice that they are symmetric. The
two infected classes have similar structures (with the exception of the parameters), as a
result, certain expressions develop into symmetric expressions. We will often find this to
be the circumstance throughout the analysis of the models. It follows from the symmetry
of the equations for I and J in the model (2-1).
Looking further into the reproduction numbers, we can obtain expressions for the rate
of each type of transmission, whether it be horizontal or vertical. We will breakdown R1.
b +ip+3i(b. -p)
b.(p + al)
b P 1- P0
bl + a1 ( (2-11)
p+acb1 p+ai& b1]
Different parts of this expression have different significance. To be precise, these
parts resemble the reproduction numbers for vertical transmission (71) and horizontal
p + au
p, + a
Notice that the reproduction number of strain one, R1, is a convex combination of
R1 and 71. So now, we will impose another parameter assumption to be held true for the
vertical transmission reproduction numbers of both strains. We will assume that vertical
transmission is of relatively small significance.
R1 = <1 (2-12)
p + ac
2 -b2 < 1 (2-13)
p + a2
Theorem 184.108.40.206. The following conditions give the i/.ir:l i; of the DFE.
a.) If RI < 1 and R2 < 1 then DFE is locally .,i-n, il'.. ically stable (l.a.s.).
b.) If Ri > 1 or R2 > 1 then DFE is an unstable saddle.
a.) For the DFE to be l.a.s., all eigenvalues of the matrix (2-9) must be negative.
We need A1 < 0, A2 < 0, and A3 < 0.
R1Z < 1 4= A2 < 0 and R2 < 1 <= A3 < 0.
Finally, using the equilibrium So, we can find A < 0.
b,(1 S*) p b1S* 1- )- p bb-)
(b, b, p i
-bx + p
< 0 (using assumption (2-6))
b.) In contrast to part (a), the DFE will be unstable if one of the eigenvalues is
positive. Thus, if either R1i > 1 or 7R2 > 1, we obtain an eigenvalue which is positive and
hence the DFE is unstable. Moreover, it is a saddle since A < 0, independent of the values
for the reproduction numbers.
2.2.3 Strain 1 and Strain 2 Equilibrium
Similar to the method for finding to DFE, we can also find the equilibrium of strain
one and strain two as Si and S2 respectively. We will start by finding S2. To find this
equilibrium, we let I* = 0 and set the derivatives equal to zero. Now, the equilibrium is
expressed as the ordered triple (S* 0, J*), where S* and J* satisfy the following equations.
0 S*[b(1 S*- J*) 2J*] (2-14)
0 = J*[b2( S* J*) ( + 2) + 32S*1 (2-15)
From this system of equations, even though it is not as easy as finding the DFE, we
can solve (2-14) for S*, and then plug the expression we get into (2-15). This will allow
us to solve for J* as an expression of only parameters and reproduction numbers. To start,
we observe that we can cancel S* from the first equation and J* from the last equation.
Then (2-14) solved for S* yields the following expression.
S* 1 (bx+ 2) J (2-16)
Now, we can substitute (2-16) into (2-15), in order to solve it in terms of P*.
0 b2(1 S* J*) (- + a2) +/2S*
/i+a2 = b2 ( +k, +032J +/2 1 b- _J b2J*
b, b, b, b )
S b( + a2)(Z2 1)
/32(b +32 2 b2)
In knowing the value for P*, we can now plug this value back into (2-16) to obtain
the value for S*. Since we have found S* and J*, this also means we have found S2, an
S( (b +32)( +2)(Z2 1) b (p + 02(2- 1)
2 (,S0J*) b, 32 /3(b, +2 b2) 2(b, 2 2 ) 18)
We can use this same process, or the symmetry of the equations, to find the
equilibrium of strain one.
-( p (b. + 13)(p + ai)(7i- 1) b.(p + a i)(7i 1)
S1 (S*,I*,0) -- bb1 b ) ,(0+b) ) (2-19)
b, b1 (b, + 01 b1) Pl(b, + 01 bi) '
Theorem 220.127.116.11. Assume R 1 < 1 and 12 < 1.
a.) S1 exists iff R7 > 1.
b.) S2 exists iff 72 > 1.
a.) In the process of finding these equilibria, we used substitution methods. For these
equilibria to exist, we need to be certain the equilibria are valid, or in other words, they
must be nonnegative. Thus, it becomes necessary for S*, J*, (1 S*), and (1 S* J*),
all to be positive.
First, it is necessary to show J* > 0. From (2-17), J* is expressed as the following.
J b(p + a2)( Z2)
2 (b2 b 2)
b.(p + 2)(72 -1) (
32(32 + bk b2)
Now we will examine the individual signs of the numerator and denominator.
b(p + a2)( (2 1) > 0 (as the parameters are all positive)
32(32 + bk b2) > 0 (by assumption (2-4))
Therefore, the overall value of J* > 0, and the entry in the ordered triple is valid.
Next, from (2-16) we have S* expressed as the following.
b, +/ 2 J
1 p b + 32 b (p + )(Z2 -1)
bx bx 02 (02 + bx b2)
1 p (b, + 02)(P + Q2) b2 (bx )
b 032(/2 + b 2) b (p + 2)
b1 2 p 2 + b b2) b (1 b}
Notice the denominator of this term is positive. Thus, in the next step, after taking
a common denominator, it is left out as it has no influence to the overall sign of the
expression. In the numerator we are left with the following expression.
S0232 + b, b2) 1 (b, + 02) !b2 +2 (lt -( +a2)1
-O22 1 (b, + 02)b2 + (b, + 02)} + Q2)
-b2 (0 + ) + (b, + 02)(/ + a2)
-b22 p(b2 b) + (P2 + a2) + b.a2
Using (2-4), we get the expression, p(b2 b,) < 0. In addition to this, since all of
the values for the parameters are positive, ba2 > 0. With these individual expressions
positive, just as we did with the denominator, we will leave them out of the remaining
32(-2 + (P + -2))
= 22 P + 2 1
S/32( + 2)(l 72)
This follows from assumption (2-13), which then implies that S* > 0.
The final expressions to check are (1 S*) and (1 S* J*), which by (2-16) are
expressed as the following.
1 S* J*
b, + 02)
+ b+ 0
+ b J
Again, we have all positive parameter values, as well as (from an earlier part of this
proof) J* > 0, which imply that (1 S*) > 0. Additionally, from this analysis, we also
observe that (1 S* J*) > 0.
Therefore, by showing these four expressions S*, J*, (1 S*), and (1 S* J*),
all to be positive with the assumption 72 > 1 and 72 < 1, we have a valid strain two
equilibrium. Hence, S2 exists.
b.) It follows analogously from part (a) to prove S, exists when Ri > 1 and 71 < 1.
2.2.4 Invasion Numbers
As mentioned before, in addition to reproduction numbers, the invasion numbers also
p11 a large role in the analysis of a particular model. Finding the invasion numbers is
similar to the process for finding reproduction numbers, with respect to using the jacobian
matrix of the system. However, we do not evaluate the jacobian at the DFE. To find the
invasion number of strain one, R we evaluate the jacobian at the equilibrium of strain
two, 2 = (S*, 0, *). We will use the eigenvalue in the column of the jacobian which was
differentiated with respect to I. The jacobian matrix (2-9) evaluated at the ordered triple
(S*, 0, *), yields the following matrix.
b( (1- 2S* J*) p/ 2J* -bJS* 3S* bSS* 32S*
0 bi(1 S* J) (p + ai) + iLS* 0 (2-22)
-b2 J + 32J* -b2 J* b(1 S* 2J*) (/ + a2) + 32S*
Notice in the second row of matrix (2-22), since we have all zero entries except
for the (2,2) entry, we can simplify the matrix around this term. Moreover, this entry
actually becomes the invasion number 7. Just as we did with the reproduction numbers,
we know this eigenvalue must be negative in order for the equilibrium S2 to be locally
.,-i~,111' I ically stable. This eigenvalue is negative if and only if the invasion number, as an
expression, is less than one.
bi(- S*-J*)- (p+ai)+A3S* < 0
bi(- S*-J*) + S* < (p + a1)
bi(1 S* J*) + QS*
p + a,
bi(1 S* J*) + 1S*
R = (2-23)
We can compute the invasion reproduction number of strain two in the same manner.
Or using that the two infected classes are symmetrical, we can determine the invasion
number of strain two by simply switching the subscripts of the corresponding parameters.
b2(1 S* I*) + 2S* (224)
p/ + a2
The biological interpretation of the invasion reproduction numbers is as follows. The
invasion reproduction number of strain one is the number of secondary cases that one
infected individual will produce in a population where strain two is at equilibrium. 71
measures the invasion capabilities of strain one. Similarly, i2 is this measure for strain
Now that we have completely determined the expressions for both invasion numbers,
we are finally in a position where we can start analyzing the model from an analytical
Theorem 18.104.22.168. Assume K2 > 1 and 72 < 1.
a.) S2 is 1, .,1,i .. 'lqil *,:. i11l' stable if 7I < 1.
b.) S2 is unstable if RI > 1.
a.) Looking at the jacobian matrix evaluated at (S* 0, J*), we get the same matrix
we used to find the invasion numbers, matrix (2-22). With zero entries in positions (2, 1)
and (2, 3), finding the determinant J AI| (where I represents the identity matrix)
becomes simpler. We can use cofactor expansion to evaluate the determinant. With this
method, the term we are expanding about, or the entry (2, 2), becomes positive. Hence,
the remaining two eigenvalues are eigenvalues of the following 2 x 2 matrix.
b(1 S*) b1S* bJ 02J* P -b-S* 02S*
-b2J* + /2J* b2(1 S*) 2b2J + S* (p + c2)
However, from the fact that S* and J* satisfy equations (2-14) and (2-15), we know
the expressions b,(1 S*) --b, 2 P p = 0 and b2(1 S*) -b2J* S* + a2) = 0.
Thus, after substitution, we are left with the remaining matrix.
-bS* -bS* 2S* (226)
-b2J* +2J* -b2J*
Det = b S*b2J* + (-b2J* + /2J*)(b S* + 02S*)
Sb1S*b2J* b2J*bS* b2J*32S* + /2J*b1S* +/2 J32S*
-b2J*32S* + /2J*b1S* +/2 J32S*
= *2S*(b b2) + 2J*S*
Trace = -bS* + (-b2J*) < 0
Thus, with Tr < 0 and Det > 0, the two eigenvalues of matrix (2-25), have negative
real part. Moreover, by the assumption 7R < 1, the third eigenvalue of the jacobian
matrix (2-22) -1 lil- negative. Finally, we have that equilibrium of strain two, S2, is locally
.-i',' .1 '1 ically stable.
b.) From the assumption that 7Ri > 1, entry (2, 2) of the jacobian matrix (2-22) is
positive. Whenever an eigenvalue of a matrix is positive, we have that the equilibrium is
Theorem 22.214.171.124. Assume R > 1 and R < 1.
a.) S, is ... .il;i '-. ,"l", i ./ll.i stable if 2-2 < 1.
b.) S, is unstable if 22 > 1.
The proof of (a) and (b) follow symmetrically from Theorem 126.96.36.199.
2.2.5 Expressing the Invasion Reproduction Numbers as Functions of the
The reason why the invasion numbers are so important lies in the fact that they
can be written as a function in terms of the reproduction numbers. For instance, we can
express Ri = F(RI, 72) and 72 = G(RI, 72). This then grants us the abilities to graph
the qualities R7 = 1 and 72 = 1 in the (RI, 72) plane and determine correct initial
parameter values that will confer the results we desire. Although this then leaves us the
task of expressing the invasion numbers as functions of the reproduction numbers.
In each of the two invasion numbers, we must first represent the identifiers for each
class in terms of parameters. So, we must find an expression for (1 S* J*) and S*.
However, these two expressions were already found in the linearization process of finding
the S2 equilibrium, as (2-21) and (2-16) respectively. Using (2-17) in the last step, we
obtain the following.
b l & 2 + J*) pl ( (bx +32J*)
( _p_) )+
b1p + py + (bj2 (31(bx+32))J*
p + a,
+ -V + ( _b z + (b ) )J
bip+3 i x-31p b+ ij2- (b+32)
p + a
bilp + 31b f031i b1/32 31(b + /32) J
bx(pI + a) bx(p + ai)
I b1/32 -/31(b + /2)
= / bx(p + al)
bl32 31(b, + 02) (p+ 0a2)(2 1)
p + a,1 02(2 + bx b2)
By a similar fashion, or again by using symmetry of the equations, we find the other
invasion number as a function of the reproduction numbers.
b20 1 32(b, +013)(p + a)(R 1)
/v2 74'2 +~~---------,-;---
R/ + ia2 /(1 1+ b bl)
In order to complete the process of expressing R71 and 7'2 in terms of R1 and 72,
we must eliminate one parameter from each of the invasion numbers entirely. These
expressions are part of R1 and 72. So we will use 7R and 72 to replace ac and a2
respectively. We can do this by solving the reproduction numbers R1 and 72 for these
bl(1 S*) + P1S*
p+ai bj=Z -
b2(1 S*) +32S*
/i+a2 b--- -
We can then substitute these expressions into the corresponding places of the invasion
numbers. Now, after simplification, we have the invasion numbers expressed in terms of
the reproduction numbers, as well as fewer parameter values than before.
b102 0ilb + 02)Z1 b2p + 32(b, I) 72 -1
'1 = R+
Sbl i + 01(b p) R2 022 +b b, b2)
SR[b 02- 1(b, + 02)][b2 +2(b-p)] 1 (2 27
=+ b1(+b 1 (2 27)
[blp + 0l(b, p)] [2 (02 + b. b2)] R2
b201 32(b + 01)R2 bi + pl(b p) R1 1
b2 +02(b, p+) R1 /1(01 + b, bl)
SR2 + 72 [b2 1- /2(b, + / 3)][bil + /1i(b p)] ( (2 28)
[b2P +/32(b, -p)] [/1(01 +b, bi)] R
2.2.6 Parametric Plot
With the help of Mathematica, a mathematical software, we can now generate
parametric plots from the invasion numbers. Having the invasion numbers in terms of the
reproduction numbers, this program provides us the capability to plot two separate
functions on the same plot. We consider the (RI, 72) plane. Since 7R1 and 7'2 are
functions of Ri and '2, the equations R, = 1 and '2 = 1 define two curves in the
(RI, 72) plane. From the first invasion number R1 = 1, we can express the reproduction
number R7 as a function of R2. Similarly, from R2 = 1, we can express R'2 as a function
To generate these parametric plots we must first set the invasion numbers equal
to one and then express them as functions of the reproduction numbers. To make this
process easier, we will take all the parameters of the invasion numbers (2-27) and (2-28),
and symbolize them as single constants ICi and /C2 respectively.
[biP2- 0l(b, + 02)[b2 + 2(b, )]
[bp3 + P0(b, )][2 2 + b b2)]
[b20 (b, + 01)][bip + A (b,- p)]
[b2P + 32(b, p)][1 ( + b bi)]
Using (2-27) and (2-28), we express 1 = 1 as R, = g(R2), and 72 = 1 as
"R2 f(R-i). Using KIC we can write the invasion number R 1 as a function of the
RI RI + Rz/Ci 1
(1 ) 1
Solving the above equation for 7 we get the consequent expression.
1 + =11 g(R 2) (2-29)
By the same process, we can express the second invasion number 7~2 as a function of
the reproduction numbers.
1 + 2(1 f(-)
Using these two functions involving the reproduction numbers, RI g(-R2) and R2 2
f((RI), we can now plot the curves in the (R1, 72) plane, which yields the parametric
plot, (Figure 2-2). This will allow us to view these functions and enable us to draw certain
conclusions. Most importantly, we will be able to determine the correct values for the
reproduction numbers which lead to specific outcomes of the system.
Before viewing the parametric plot, we will analyze the two curves. If we can
guarantee certain situations arise within the parametric plot, this will then correspond
to the correct conditions needed to show the dominance of a strain with suboptimal
reproduction number. To show that this situation will occur, we will examine the constant
values of ICI and /C2, which encompass all the parameter values of the invasion numbers.
From Theorem 188.8.131.52, we know that if 72 < 1, then strain two cannot invade the
equilibrium of strain one. In other words, strain one dominates. Moveover, strain one
has the ability to dominate with a smaller reproduction number if R2 > R1, or when
f(7(R) > 7R. In order to get this situation on a parametric plot, we need the curve f7(Ri)
to be above the line R1 = 72. For this to occur we need to satisfy the following condition.
(1 > 7 (2-30)
Notice that f(7R) is a decreasing function of R7 if /C2 > 0, and f(7R) is an increasing
function of i, if C2 < 0. Also, if C2 > 0, the largest value of f(Ri) is f(1) 1 =
f(RI) < Ri for all R1 > 1. So the condition in (2-30) will never hold.
Now, take C2 < 0 = f (R1) is an increasing function of 1i. So f(1) 1 = f (R1) > 1
for all R1, so the condition in (2-30) may hold. There are two distinct cases we will
examine to determine exactly when (2-30) holds.
In the first case, we will assume 1 + C2 > 0. Then 1 + K2(1 >) > 0 for all R1 > 1.
To see this, notice that 1 + 2(1 -)2 is a decreasing function of 71, so the smallest value
is when R1 i oo. To find this value explicitly, we evaluate the following limit.
lim 1 + C2 2 ) 1+/C2
(- (oo ( 11R ;
If 1 + 2Ca > 0, then f(R ) is a continuous positive function of R1 for all R1 > 1. From
this, we can see directly that (2-30) will not hold.
1+ (1 ) > 2
1 > 1 + 2 (1 ]
SR1 + C2(1- 1)
S(1 +C2)71 1C2
This is a convex combination of R1 and one, so its value is alv-iv- greater than one
and we will never be able to attain the condition.
In the second case, we assume 1+/C2 < 0. But, with this assumption, the denominator
of f(7R) may become zero. This is because for some value of 7R, denoted by R7, the
equality 1 + C2(1 ) = 0 holds. This value can be calculated in the following steps.
S/ 12 +1
Clearly, R7 > 0, as both the numerator and the denominator are negative. Moreover,
R* > 1, since the numerator is larger than the denominator, which holds in the following.
lim > n\
Rl-. 1 + / t-(1 -2 I
Therefore, there are values of R1 for which f(7R) > 7R. Thus, the inequality (2-30)
It is understood, as we increase the values of R7 and 72 in the parametric plot, the
value of one of the curves will switch from positive to negative. By evaluating the above
limit, we see that (2-30) will hold, but only for some finite R1, specifically when R7 < R4
(or the value where the plot switches from positive values to negative values).
With a better understanding of the parameter situations needed to generate a specific
plot, we will now examine the parametric plot in (Figure 2-2). To obtain the specified
parameter values, we generated a random search until all given assumptions were met.
12 Parametric Plot
3 1=1 /.,= 1
0.5 1 1.5 2 2.5 3
Figure 2-2. Parametric plot of the invasion numbers R1 and R2 expressed as the
reproduction numbers. Parameter values used for this plot are, b, = 0.7791,
bi 0.4261, b2 = 0.6864, /1 = 0.2396, /2 0.765, and p = 0.0713. These
parameter values yield Ci = -0.19 and /C -1.47.
To distinguish the importance between each region of this plot, we will shade various
sections differently in (Figure 2-3). The white area of this plot is not useful for us. The
reason for this comes from the fact that in this area, both reproduction numbers are less
than one. This then implies, in most cases, that both strains will die out.
One of the more important areas of this plot is the darker gray area (in the upper left
of Figure 2-3). Here we have the invasion number R2 > 1, while R1 < 1, so strain two has
2 Fully Shaded Parametric Plot
0.5 1 1.5 2 2.5 3
Figure 2-3. Shaded parametric plot of Model 1 to distinguish between different areas.
the ability to invade strain one. For this reason, reproduction numbers with values within
this area will yield a system in which strain two will be the dominant strain.
Likewise, in the light gray area (towards the right of Figure 2-3), we have a similar
scenario with the other strain. In this area 72 < 1 and 1 > 1, which leads us to conclude
strain one will dominate.
The third area of importance is the black region. Here we have that both invasion
numbers are greater than one, or that either strain has the ability to invade the other
strain. However, since both invasion reproduction numbers are greater than one, we expect
that coexistence of the two strains will occur. So for the reproduction numbers, R7 and
72, whose values fall within this area, both strains will coexist within the population.
With all of the different shades explained, we now must focus on the most important
information that this parametric plot offers us. In addition to plotting the two functions
of the invasion numbers, we have also plotted a few other lines to serve as guides. In
particular, the diagonal line gives the value where R1 = 72. In any case of dominance,
we should be able to conclude that the strain with the higher reproduction number
will dominate the strain with the lower reproduction number. This is the rational in
determining which of the strains dominates in a particular section. But, if we look closer
at the light gray area, we notice something very important.
12 Shaded Parametric Plot
1 - - - - - -
0.5 1 1.5 2 2.5 3
Figure 2-4. Parametric plot highlighting the area of greatest importance.
(Figure 2-4) isolates a very important region of the original parametric plot (Figure
2-3). As mentioned before, within the entire light gray region of (Figure 2-3), strain one
dominates. But, in the light gray area of (Figure 2-4), we also have R1 < K 2. Thus, we
have found a region where a strain with suboptimal reproduction number can dominate.
Our hypothesis states that in this case, when vaccination is introduced, strain replacement
To demonstrate the strain replacement effect, we will produce a plot, (Figure
2-5), of system (2-1). We will choose the parameter values in such a way that the two
reproduction numbers will fall within this shaded area. Thus, to show strain replacement,
we will use a point labeled "D" corresponding to an ordered pair (RI, 72) within the
In (Figure 2-5) we have the reproduction number of strain one, R, = 2.00, while the
reproduction number of strain two, 72 = 2.51. Yet, strain one, represented by I(t), is
the dominant strain, while strain two, represented by J(t), is eliminated. Therefore, the
Suboptimal Reprodcution Number
20 40 60 80 100t
Figure 2-5. Simulation of Model 1 when the strain with lower reproduction number, I(t),
dominates a strain with a higher reproduction number, J(t). Parameter values
used for this plot are, b, = 0.7791, b1 = 0.4261, b2 = 0.'. 1 = 0.2396,
32 0.765, p = 0.0713, ac = 0.0178, and a2 0.1382. The reproduction
numbers follow as, R = 2.00 and Re2 2.51.
principal that the strain with higher reproduction number excludes the strain with lower
reproduction number is no longer valid in this model .
Using the parametric plot, (Figure 2-2), we were able to show a particular situation,
namely dominance of a strain with suboptimal reproduction number, can be attained. But
just as importantly, we can also use the analysis to choose parameter values which do not
provide us with a certain situation. To generate these plots, it must be the case that one
of the assumptions on the parameter values, (2-3), (2-4), and (2-6) is no longer held. For
this example, we will not assume (2-3), or that b, / b1. With this in mind, we can find
CI = 1.75 and C2 = 5.62, so both values are positive and generate the plot, (Figure 2-6).
Unlike (Figure 2-2), in (Figure 2-6) there is no region which sil-.-. -I- that a strain
with suboptimal reproduction number dominates a strain with a higher reproduction
number. This is because in (Figure 2-6) the area underneath the 72 = 1 curve would be
dominated by strain one. Also, the region to the left of 71 = 1 curve would be dominated
by strain two. Finally, the region between the two curves represents coexistence. Thus,
there is no such area as in (Figure 2-4), where a strain with a lower reproduction number
would dominate a strain with a higher reproduction number.
R2 Parametric Plot
1 -- -- -- -- -- -- -- -- -- -- -- -- -- ---- -- -- -
0.5 1 1.5 2 2.5 3
Figure 2-6. Parametric plot of Model 1 when we do not get the situation of the strain with
suboptimal reproduction number being dominant. Parameter values used for
this plot are, b, = 0.423, bl = 0.5945, b2 = 0.3808, 31 = 0.0572, /3 0.1094,
and p= 0.1458
2.3 Model 2: A Model with Healthy Births from Infected Individuals
Contrary to Model 1, when we permit healthy newborn individuals from the two
infected classes I(t) and J(t), the parameters r,1 and r72 are not equal to zero. The analysis
of this model will prove to be more difficult as there are more variables to consider and
2.3.1 Equilibrium and Reproduction Numbers
Utilizing the same methods as with Model 1, we compute the DFE of this model to be
an ordered triple.
so --b 0 0) (2-31)
0 bx b
Likewise, the reproduction numbers of this model arrive to be the following.
bi + 3(b p) b2+ 2(b p) (2
R1 22 (2-32)
b(p + ari) b(p + a2)
It just so happens that the DFE and reproduction numbers are the same for Model 1
and Model 2, but these expressions will not ahv--, be of this particular form.
2.3.2 Strain 1 and Strain 2 Equilibrium
Just as we did with Model 1, we want to compute the ordered triples Si = (S*, I*, 0)
and 2 = (S*, 0, J*). Furthermore in finding S2, we set the derivatives equal to zero and
solve for the variables in the system.
O = (bS* + T2J*)(1 S* J*)- S* 02J*S*
0 0 (2-33)
o J* [b2(1- *) ( + a2) +/30S*]
Solving the last equation for either S* or J* will allow us to substitute an expression
into the first equation to solve for the other variable. We start by solving the last equation
S* -) + (/ +a2) f (*) (2-34)
Additionally, we note that if 32 > b2, then f(J*) is a linear increasing function of J*.
We compute the values of this function at zero and one as the following.
-b2 + (f/ + a2)
2 b2 (2-35)
Sf(1) / + 2 (2-36)
Due to the complexity of evaluating this expression, we will substitute f(J*) into
the first equation of the system and not the entire expression for S*. This will make
the analysis more simple when dealing with the quadratic expression we acquire after
substitution into the first equation of (2-33). In doing so, we obtain an equation for J*
and we will call it F(J*) 0.
F(J*) = (bf(J*) + /2 *)(1 f(J*) *) f(J*) P*f(J*) (2-37)
Notice that F(J*) is a continuous function of J* and the values of this function at
zero and one are as follows.
SF(0) = (bf(0))(1 f(0)) p f (0)
f (0)[b(1 f(0)) ] (2-38)
SF(1) (b f(1) + 92)( f(1)- )-f()-32f(1)
-f( 1)', f() + r2 +P +/32] (2-39)
With the right hand side of the first equation in the system (2-33) expressed as a
function of J*, we can prove that the strain two equilibrium exists. Moreover, we can show
that the equation F(J*) = 0 has zero, one, or two solutions using the Intermediate Value
Theorem. The equation F(J*) = 0 is a quadratic equation. Although we can compute J*
from it, the expression will be complicated. Without an explicit expression for J*, there
is no way for us to determine the specific equilibrium S1 or S2. However, to prove the
equilibrium actually exists, we will incorporate the following theorem.
Theorem 184.108.40.206. Assume if 72 < 1 and 72 > 1, then there is a unique solution for *.
From the assumptions .2 > 1 and 72 < 1, we get the inequality 32 > b2. This then
implies that f(0) > 0 from the expression in (2-35).
In order to use the Intermediate Value Theorem, we must test the individual signs of
F(0) and F(1).
F(0) f(0) [b (1 f(0) ]
[b (1 f(0)- ] (since f(0) > 0)
b + a2)[2 -1] > 0 (by assumption R2 > 1)
Also, from (2-38) and (2-39), we have that F(1) < 0. Hence, from the Intermediate
Value Theorem, there is a unique solution for J*.
To show this unique solution visually, we can generate a plot of the function F(J*)
from (2-37). If we make certain the assumptions from Theorem 220.127.116.11 hold, we get a plot,
(Figure 2-7), where F(0) = 0.09, F(1) = -2.08, with a zero at J* = 0.24.
F(J*) Existence of a zero
0.2 .4 0.6 0.8 1
Figure 2-7. Plot showing existence of a unique solution of J* using the Intermediate Value
Theorem. Here, it is necessary that the conditions of Theorem 18.104.22.168 hold, so
the value of 7'2 0.75 and 72 = 1.58. Parameter values used for this plot are,
b, = 0.9404, bi 0.1123, b2 = 0.3672, P1 = 0.2287, /2 0.9532, = 0.2901,
r1i 0.4355, l2 = 0.4682, a 1 0.151, and 2 = 0.198.
Now, if we assume the conditions of Theorem 22.214.171.124 do not hold, we get a quadratic
with no solution on the interval (0, 1). So assuming 72 < 1 and 32 > b2, we can show that
there may be no equilibria, as in (Figure 2-8).
2.3.3 Invasion Numbers
Establishing the invasion numbers of Model 2 follows precisely the same methodolgy
as that of Model 1. We compute 7R1 and 7'2 as the following expressions.
bil( S* J*) + 1S* b2(1 S* J) + 0S*
R =< 1 22 <1
(p + ai) (p + 2)
Moreover, we can use the expression for S* in (2-34), to substitute an expression
for S* and (1 S*). Applying this substitution, along with simplification, we obtain the
F(J*) No Solution
0.2 0.4 0.6 0.8 1 J
Figure 2-8. Plot verifying there may not be a solution for J* (using the Intermediate Value
Theorem), if the conditions from Theorem 126.96.36.199 are not satisfied. Here, the
value of 722 0.78 and R72 0.81. Parameter values used for this plot are,
b, = 0.6099, bi 0.3222, b2 = 0.5826, P31 0.5469, 02 0.6056, = 0.012,
T1i 0.7323, 92 = 0.6908, aci 0.6447, and a2 = 0.732.
invasion number of strain one as a function of J*, and of strain two as a function of I*.
S (3b2 blP2)J* + bl2 lb2 b(p+ 2)+ ( 2) (240)
(32 b2) (/i + a 1)
&2 (302b b2/1)I* + b231 32b bi(p + a) + 2( + a) (2 4)
(si bi)(p +a2)
2.3.4 Simulation Plots
Without the parameters 711 and T/2 equal to zero, the analysis for this model is difficult
to carry out mathematically. Determining initial parameter values would prove to be
a very difficult task. Even if conditions were established, due to the complexity of the
model, there is no -,iwing if the conditions were accurate. This is where we must rely
on another source, in simulations, to provide evidence of the dominance of a strain with
suboptimal reproduction number.
Through simulations of the model, we can produce different results of the system
with given parameter values. For instance, when we choose parameters which yield
reproduction numbers of the two strains I and J as, R, = 2.53 and 72 = 3.79 respectively,
we can then generate a simulation showing how strain J dominates. This is due to the fact
that the strain with the higher reproduction number should dominate a strain with a lower
reproduction number. This is shown through the first simulation, (Figure 2-9).
Simulation Model 2
10 20 30 40 50 60 70
Figure 2-9. Simulation of Model 2 when the strain with a larger reproduction number
dominates. Parameter values used for this plot are, bx = 0.8618, bl = 0.8031,
b2 = 0.6802, Pl = 0.8833, 32 = 0.9589, = 0.1144, r = 0.188, 72 = 0.2847,
ac = 0.2294, and a2 = 0.1286.
In contrast to (Figure 2-9), we can also simulate the model to generate a different set
of results that we are interested in. With the second simulation, (Figure 2-10), we have
chosen parameters that yield the reproduction numbers of strains I and J to be R = 2.88
and R2 = 3.62 respectively. Conflicting with the results from the first simulation, we now
have shown the strain with suboptimal reproduction number, strain I, dominates strain J,
which has the higher reproduction number.
Suboptimal Reproduction Number
0.2 '.. S(t)
50 100 150 200 250 300 350
Figure 2-10. Simulation of Model 2 when a strain with suboptimal reproduction number
dominates. Parameter values used for this plot are, b, = 0.7791, b1 = 0.4261,
b2 = 0.6864, P31 0.2396, /2 0.756, p = 0.0713, rl = 0.1923, 92 = 0.6469,
ac = 0.0178, and a2 = 0.1382.
VERTICAL TRANSMISSION MODEL WITH VACCINATION
In this chapter we will examine another model, one which is responsible for answering
the main question of this study. This model is very similar to the model in chapter two,
however, a vaccination class is now included.
3.1 Model with Vaccination
Below is the general model that includes a vaccinated class. The vaccination within
this model is assumed to be "perfect", or that it is one-hundred percent effective against
both strains in the system. For the meanings of the parameters, refer to (Table 2-1).
d= (b1S + 711 + T2J + r/3V)(1 S I J V) (p + )S P11S 02JS
d =I[b(1 S- I J- V) (+ ai) + 3S]
= J[b2(1 S I J V) ( + a2) + 2S]
dt S V
The flow diagram for this system is the following.
Figure 3-1. Flow chart of the model with vaccination.
The same assumptions (2-3), (2-4), and (2-6), which were used for the model in
chapter two, will also be assumed for this model. We will make one additional assumption.
b, = /3 (3-2)
This assumption signifies that the per capital birthrate of susceptible newborns into
the vaccination class is equal to that of the susceptible class.
3.2 Model 3: A Vaccination Model without Healthy Births from Both
For this model, T11 and 92 are assumed to be zero. The reproduction numbers of the
two strains are the following.
bi(-) + 31i(77)(1 ) b2(') + 02(.')(1 )
bi() P b bx P- 72b b3)
Ip + ac i p + ba2
With a vaccinated class, the reproduction number values depend on the value of Q. In
fact, the reproduction numbers are decreasing functions with respect to Q. In the absence
of vaccination, = 0, the reproduction numbers are identical to the reproduction numbers
(2-10) in chapter two. The system above albv--- has a disease-free equilibrium.
s ((t b. + b. + \
Through computation involving the linearization process, we are able to find values
for I* and *.
P (hlt 3 1 +)[ 1 PR ()J 1[1 l (R) 1 (3 4)
3.[1 +p[_p1 C2
P+L b: (p1+) b C2
(bA1 i + P t [ ) 1 T)[1 1
P* bx P+ b +R 2)\ b R2 () (3-5)
W\[1 + p2 b2 E ) D2
P+L bz (p(+ ) bz
The expressions for I* and J* become important substitutions into the invasion
numbers, which are computed as the following.
For simplification reasons, we have taken certain groups of parameters in the
expressions for 7Z1 and 7R2 and represented them as constants. Since all of these parameter
values are constant, we can express all parameters in terms of a single constant. The
following shows the constants M1 and M2, used for simplification in R71.
bi I 3b1p !3ip I (hftbP 0p 02 *
p1 + +
SM2 ( 2 (- )
(2 + p---- (p + -) (/b 1 b. p--)- (3-7)
P + a2
For simplification reasons, we have taken certain groups of parameters in the
expressions for 71 and 72 and represented them as constants. Since all of these parameter
valuThe same cnstan be can express all the parameters in terms of a sing the constant and With
thesefollowing substitutions, the ionsta numbers can be expressed more simplification in .
1 = 1 2J* (3 8)
p + a /
b p bp + +2
M2 p (b, A/0 2
(p+ b. b.1 bp+ )
p + a2
-&2 -V --------
It is now necessary to generate the parametric plots for this model. To do so, we
must find the expressions 7Zi(n) = F(71i(), 72()) and 72() = G(Ri(), R2 ())
After we have these expressions, just as we did in chapter two, we will need to represent
the reproduction numbers as functions of one another. For example, we will show how to
obtain 7R,() ( f(7R2(2)), while the second expression will follow by symmetry.
To start, we need to set the invasion number 7i = 1 and then find the necessary
expression. After substituting the values for J* and (p + ai) from (3-5) and (3-3)
respectively, into (3-8), we obtain the following expression.
I (4, + 422)
W) 1 Q 1____
2 D2 1
M1 [ + 2M0 1()
Again, since Mi and Di are constants (with i = 1 or 2), we can simplify the
expression even more. To do so, we define a new constant Ci, where IC = M1 and
K2 = 2DIM Thus, the expression RZ1(Q) f (2(S)) is as follows.
R () K + 2 1 (3-9)
ICI + /2[i R 2yb
Similarly, this same process is used to find Re2() g= (Ri()) with the exception that
the constants will now be labeled as Li.
R2 ( ) L (3-10)
+i 2[l l+ 2y)
Now, using these two functions, we can generate a parametric plot of the two
curves. However, since we now have vaccination in the model, these plots will fluctuate
depending on the level of vaccination. We will use this to demonstrate how specific areas
of dominance within the parametric plot change with increasing levels of vaccination.
Simulation Results. The following simulations of Model 3 provide evidence of
strain replacement with increasing levels of vaccination. With each of the simulation plots,
there exists a corresponding parametric plot. However, since the reproduction numbers
are decreasing functions of i, increasing the level of vaccination changes these values. In
addition to this, the parametric plots will also change as they utilize the values of the
Given a set of parameter values, we can use (3-9) and (3-10) to generate a parametric
plot. From the given parameters, along with the level of vaccination, we can compute
the specific reproduction numbers Ri and 72. With these values, we then obtain an
ordered pair (RI, R2). Depending on the location of this point within the parametric plot,
K2 Parametric Plot with =0
RI=1 Model 3 =O0
6 ." 0.5 I(t)
2 : /0' 0 2 ". .
.1 ".: ..- S (t)
1 2 3 4 6 7 1 20 40 60 80 100 120 140
Figure 3-2. Parametric and corresponding simulation plot of Model 3 when =- 0. The
reproduction numbers of I and J are, 71 = 4.36 and 72 = 7.05 respectively.
(a) Parametric plot. (b) Simulation.
will give rise to particular situations of the system. These reproduction numbers are the
coordinates of the point "D" in the parametric plot.
The parameters used for this simulation are, bx = 0.732, b = 0.197, b2 = 0.0652,
31 = 0.440, /32 0.876, p = 0.0771, Th = 2 = 0, 3 = b,, ac = 0.0179, and a2 = 0.0349.
In the first of three series of simulation plots, (Figure 3-2), for this model, we take
vaccination to be zero, or =- 0. To understand how the system may react, we generate
the parametric plot, (Figure 3-2(a)), using the fixed parameters. Even though this
parametric plot is not an exact copy as in Section 2.2.6, it is similar and offers us the
opportunity to reference (Figure 2-3). Now, it will be easier to understand the distinct
areas where a particular strain will dominate.
The positioning of the point "D" in the parametric plot, (Figure 3-2(a)), -r-'-I -r; that
a strain with suboptimal reproduction number will dominate. We confirm this conjecture
with the outcome of the dynamic simulation, (Figure 3-2(b)), of the system.
As we increase vaccination to a level of = 0.15, not only do the reproduction
numbers change, but also the position of the point "D". Now this point falls in a region
where the strains coexist. This is shown in the second simulation, (Figure 3-3).
6 2 Model 3 #=0.15
1 2 3 4 5 6 7 0.5 1 15 2 2.5 200 400 600 800 1000
(a) (b) (c)
Figure 3-3. Parametric and corresponding simulation plot of Model 3 when = 0.15. The
reproduction numbers are, R = 1.63 and 72 = 2.44. (a) Parametric plot. (b)
Enlarged parametric plot. (c) Simulation.
R2 Parametric Plot with =0.28
1=1 2=1 2 parametric Plot with l =0.2
7 Model 3 i=0.28
4 / 0 .3 .
3 : 0.75 .0 .2 "
2 : / 0.5
1 0. 25 0.1 S(t)
S 2 3 4 5 6 l 0.25 .5 .75 1 188.8.131.52 200 400 600 800 1000
(a) (b) (c)
Figure 3-4. Parametric and corresponding simulation plot of Model 3 when = 0.28. The
reproduction numbers are, R = 1.11 and 72 = 1.57. (a) Parametric plot. (b)
Enlarged parametric plot. (c) Simulation.
Finally, in the third simulation plot, (Figure 3-4), we increase vaccination even further
to a level of = 0.28. Even though there is not much change in the parametric plot,
(Figure 3-4(b)), the location of the reproduction number coordinate "D" proves to be very
important. Now, this point falls into a region where the opposing strain, one which was
not dominant with zero vaccination, will now be the dominant strain.
Through this simulation we see that strain replacement occurs in a vaccination model
where 1, i11,!: births from both infected individuals are not permitted.
3.3 Model 4: A Vaccination Model with Healthy Births from Both Infected
For this model, we are utilizing the entire system (3-1).
Simulation Results. The following simulation demonstrates how strain replacement
does in fact occur with perfect vaccination through vertical transmission. Although it is
possible to generate the parametric plots for this model, the computation would prove to
be very difficult. If the parametric plots were generated, we would see the same result that
we obtained in Model 3. Increasing the level of vaccination would change the positioning of
the reproduction number coordinate "D", thus leading to certain simulation scenarios and
The parameters used for this simulation are, b, = 0.7791, bl = 0.4261, b2 = 0.6864,
3i 0.2396, 32 = 0.765, p = 0.0713, 1 = 0.1923, 92 0.6469, 3 = b,, a, = 0.0178, and
2 = 0.1382. Only the value of i will vary.
The first simulation plot, (Figure 3-5), has a vaccination level of zero, ) = 0.
Through this simulation we observe that a strain with suboptimal reproduction number, I,
dominates strain J, which has a higher reproduction number.
Model 4 = 0
0.2 \ S(t)
100 200 300 400 500
Figure 3-5. Simulation of Model 4 when = 0, showing how strain I is dominant even
though it has a lower reproduction number. The reproduction numbers of I
and J are, Ri = 2.88 and R2 = 3.61 respectively.
As vaccination is applied to the population and raised to a level of = 0.2, we have
coexistence of the two strains. This is seen through the simulation plot, (Figure 3-6).
Model 4 *=0.2
0.1 --- ---- --
J t I(t)
500 1000 1500 20002500 3000 3500
Figure 3-6. Simulation of Model 4 when p = 0.2, showing coexistence of the two strains.
The reproduction numbers are, R = 1.07 and 72 1.17.
Finally, as vaccination is increased more, to a level of p = 0.23, we obtain the
simulation, (Figure 3-7), where strain J dominates strain I. But recall, in the first
simulation, where p = 0, we had that strain I was the dominant strain. This is where
strain replacement occurs.
Model 4 i=0.23
0 .i \ _I (t )_ J t
1000 2000 3000 4000 5000
Figure 3-7. Simulation of Model 4 when p = .23, showing the occurrence of strain
replacement. The reproduction numbers are, R1 = 1.02 and 72 = 1.08.
The final simulation, (Figure 3-7), verifies that strain replacement occurs in a
vaccination model where newborns from both infected individuals are allowed. Also,
since the reproduction numbers are decreasing functions of i, as we increase the level of
vaccination within the population, the reproduction numbers decrease. Eventually as we
increase the level of vaccination enough, both strains will be eradicated.
The main topic of this study was to investigate a certain hypothesis in regards to
perfect vaccination as a mechanism for serotype (strain) replacement. It has already been
established, in other mathematical papers , , and , that strain replacement
occurs when the vaccine is differentially effective, as well as in specific cases when it is
perfect. If the vaccine is perfect, for strain replacement to occur, a trade-off mechanism
should also be in place. Below is a table which summarizes all of the different trade-off
mechanisms and if strain replacement with perfect vaccination occurs with the given
Table 4-1. Trade-off mechanism overview
Strain with lower repro- Strain replacement with
Trade-off mecha- duction number excludes perfect vaccination
nism strain with higher repro-
Model 3 Yes Yes
Model 4 Yes Yes
Super-infection Yes Yes
Co-infection Yes Yes
Mutation No No
Cross-immunity No No
The analysis and simulations of the four models in this study provide us the evidence
needed to prove the hypothesis dealing with perfect vaccination. The results from Model
3 and Model 4 are of particular importance to us as they deal with vaccination. However,
since Model 1 and Model 2 function similarly to the model with a vaccinated class, we
were able to gain understanding on parameter values and other assumptions needed to
achieve the desired results.
The hypothesis presented sl:.;; -1 that only trade-off mechanisms which lead to the
dominance of a strain with suboptimal reproduction number lead to strain replacement
with perfect vaccination. This hypothesis is confirmed to be true using the simulations
involving the model with vaccination. In both cases of the vaccination model, with and
without healthy births from infected individuals, we can conclude that perfect vaccination
leads to serotype replacement.
To address this hypothesis further, we must show that when the vaccination level
is zero within the vaccination model, it must be the case that a strain with suboptimal
reproduction number can dominate within the population. This is where the analysis from
Model 1 and Model 2 are important to us. When there is zero vaccination, the systems
in reference to these models represent the vaccination model with zero vaccination. Thus,
as proved in Section 2.2.6, there are values of the reproduction numbers which permit a
strain with suboptimal reproduction number to dominate.
After proving that a strain with suboptimal reproduction number can dominate
with zero vaccination, it is then somewhat intuitively clear why this is an essential
component to show strain replacement occurs. Without this criterion, the strain with
higher reproduction will dominate with zero vaccination and then it will continue to
dominate regardless of the vaccination level.
On the other hand, with the capability for a strain with lower reproduction number
to dominate when the level of vaccination is zero, we are able to get a replacement effect.
Since both reproduction numbers are decreasing functions of the vaccination level i, as
the level of vaccination increases, the reproduction numbers of the two strain decrease.
As this happens, not only does the parametric plot change, but also the location of the
reproduction number coordinate point on the plot.
We witnessed with the simulation of Model 3, that the coordinate signifying the
value of the reproduction numbers falls into different regions of the parametric plot as
vaccination is increased. Moreover, as it changes location it produces different dynamic
simulations of the system. This point starts in a location where a strain with suboptimal
reproduction number dominates. Then at a certain level of vaccination, the point -.: -
that the two strains are able to coexist. Finally, increasing vaccination levels further, leads
to strain replacement. This is because the dominant strain with zero vaccination dies out,
while the other strain persists and becomes dominant.
This same result follows for Model 4, where healthy births from both infected
individuals are permitted. Hence, if a strain with suboptimal reproduction number
has the ability to dominate with zero vaccination, then vertical transmission as a trade-off
mechanism, will permit strain replacement through perfect vaccination.
 A. M. Behbehani, The smallpox story: Life and death of an old disease,
Microbiological Reviews, 47 (1983), 455-509.
 H. Bremermann and R. H. Thieme, A competitive exclusion principal for pathogen
virulence, Journal of Math. Biol., 27 (1989), 179-190.
 S. Busenberg and K. Cooke, Vertically Transmitted Diseases: Models and Dynamics,
Biomathematics 23, Springer, Berlin, 1993.
 J. W. Dale, Molecular Genetics of Bacteria, John Wiley and Sons Ltd., Chinchester,
 A. d'Onofrio, G'1,/l.l,.,1 stable vaccine-induced eradication of ho,.:...,.,Ill, and verti-
S'll/; transmitted infectious diseases with periodic contact rates and disease-dependent
demographic factors in the population, Applied Mathematics and Computation, 140
 M. El-Doma, Aal;,-.:- of an SIRS age-structured epidemic model with vaccination
and vertical transmission of disease, Applications and Applied Mathematics Journal,
1 (2006), 36-61.
 D. A. Henderson, The looming threat of bioterrorism, American Association for the
Advancement of Science, 283 (1999), 1279-1280.
 D. R. Hopkins, Smallpox: Ten ;;, r,- gone, American Journal of Public Health, 78
 T. Horimoto and Y. Kawaoka, Pandemic threat posed by avian :,fl. ,. ../ A viruses,
Clinical Microbiological Reviews, 14 (2001), 129-140.
 M. Iannelli, M. Martcheva, and Xue-Zhi Li, Strain replacement in an epidemic model
with super-infection and perfect vaccination, Math. Biosci., 195 (2005), 23-46.
 H. Inaba, Mathematical ,,.ir ..:- of an age-structured SIR epidemic model with
vertical transmission., Journal of Discrete and Continous Dynamical Systems (Series
B), 6 (2006), 69-96.
 J.-H. James Ou, Hepatitis Viruses, Kluwer Academic Publisher, Massachusetts, 2002.
 S. A. Kauffman, The Origins of Order, Oxford University Press, New York, 1993.
 C. M. Kribs-Zaleta and M. Martcheva, Vaccination strategies and backward
bifurcation in an age-since-infection structured model, Math. Biosci., 177-178
 M. Lipsitch, S. Siller, and M. A. Nowak, The evolution of virulence in pathogens with
vertical and horizontal transmission, Evolution, 50 (1996), 1729-1741.
 M. Martcheva, On the mechanism of strain replacement in epidemic models with
vaccination, Current Developments in Mathematical Biology Proceedings of the
Conference on Mathematical Biology and Dynamical Systems, World Scientific
Tyler, The University of Texas at Tyler, 2007.
 M. Martcheva, B. M. Bolker, and R. D. Holt, Vaccine-induced pathogen 1,''
replacement: What are the mechanisms?, submitted.
 R. May and M. Nowak, Coinfection and the evolution of parasite virulence, Proc.
Royal Soc. London B, 261 (1995), 209-215.
 M. Nowak and R. May, Superinfection and the evolution of parasite virulence, Proc.
Royal Soc. London B, 255 (1994), 81-89.
 D. E. Swayne and others, C'l,/.,. I, rization of highly pathogenic H5N1 avian
:,lI;,. ..i A viruses isolated from South Korea, Journal of Virology, 79 (2005),
 H. R. Thieme, Mathematics in Population Biology, Princeton University Press,
Princeton and Oxford, 2003.
Douglas Harold Thomasey was born on June 25, 1982 in Somerville, New Jersey. He
spent the first eighteen years of his life in the small town of Middlesex with his parents
Richard and Laura, as well as his two older brothers Bryan and Kevin. Family life was
very important to Doug, and he would consider his family to be very close and loving.
After graduating from Middlesex High School, Doug completed his undergraduate
studies at Lynchburg College in Virginia. While at Lynchburg, Doug was a varsity athlete
in track and field specializing in the decathlon. With his extremely hard work ethic and
determination, Doug achieved All-American status in 2004 and 2005 by placing 5th and 6th
in the nation respectively. It was also on the track where Doug met his girlfriend Ashley
Palmer, who happens to be a 3-time All-American.
In addition to the demands of a grueling training regiment in track and field, Doug
also found success academically, double in, Pi ii.w'-; in both mathematics and Spanish.
The results of his thesis in mathematics entitled, "The Effects of a Parasitic Copepod
(achtheres) on Smith Mountain Lake", was presented before the Virginia Department
of Game and Inland Fisheries, and also received honors at Lynchburg College's Annual
Student Scholar Showcase.
With these combined successes in academics and athletics, Doug won the Lynchburg
College Scholar-Athlete of the Year award for the years 2003 and 2005. More notably,
he received a NCAA Postgraduate Scholarship, given to students who excelled both
academically and athletically at a national level. He graduated in May 2005, Magna Cum
Laude with a 3.87 grade point average.
With some influence from his undergraduate coordinator, Dr. Kevin Peterson,
an alumni of the University of Florida, Doug decided to attended graduate school for
mathematics at UF. After his first year of study, he realized that his greatest interests
lie within the field of biomathematics. While in a biomathematics course his first year,
Dr. Martcheva gave a talk that interested him. He then approached her in regards to a
possible topic for his thesis and choose her as his advisor.
In completing his M.S. degree at the University of Florida in May 2007, Doug hopes
to secure a job in teaching. He will further explore his options in the future with respect
to completing a Ph.D., as well as different career paths, but for the present he hopes to
pursue his athletic goals, along with Ashley, in track and field.