A N EPIDEMIC M O D E L W I T H I M M I G R A T I O N OF INFECTIVES A N D V A C C I N A T I O N By Eunha Shim B.Sc. (Mathematics) The University of British Columbia, 2002 A THESIS S U B M I T T E D IN P A R T I A L F U L F I L L M E N T O F T H E R E Q U I R E M E N T S F O R T H E D E G R E E O F M A S T E R O F S C I E N C E in T H E F A C U L T Y O F G R A D U A T E S T U D I E S D E P A R T M E N T O F M A T H E M A T I C S I N S T I T U T E O F A P P L I E D M A T H E M A T I C S We accept this thesis as conforming ; to the required standard THE UNIVERSITY OF BRITISH COLUMBIA APRIL 2004 © Copyright by Eunha Shim, 2004 Library Authorization In presenting this thesis in partial fulfillment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Shim, Eunha 20/04/2004 Name of Author (please print) Date (dd/mm/yyyy) Title of Thesis: An epidemic model with immigration of infectives and vaccination. Degree: M.Sc. • Year: 2004 Department of Mathematics The University of British Columbia Vancouver, BC Canada Abstract This paper is an investigation of the possible behaviors of epidemics when there is immi-gration with some infectives and when vaccination is in effect. Infective immigrants from outside the population may introduce communicable disease into a host population. To help controlling disease one might introduce mass vaccination. First the model of SVIS type disease is discussed to describe the behavior of an epidemic disease when spreading into a population with immigrants and when a vaccination policy is in effect. Even in this simple version of the model, backward bifurcation and multiple endemic steady states can be observed with some sets of parameter values; by mathematical analysis the condition for a stability change and for having a backward bifurcation is proved to be identical in this case. Also the sufficient condition for backward bifurcation is stated explicitly. In the case of forward bifurcation, there can be exchange of stability and Hopf bifurcation occurs in such cases - numerical examples are provided to help understanding. Next we assume a disease with no immunity gained by recovery from infection, so the model is of SVIR type. In the case of no disease fatalities several special cases are discussed; stability analysis of steady states in each case and threshold values determining existence of endemic equilibria are presented if there are any. The proof that backward bifurcation cannot occur with any set of nonnegative parameter values in SVIR model is presented. On the other hand some numerical examples of backward bifurcation with negative immigration rates, i.e. outgo-ing immigration from host population, are shown. When there is a backward bifurcation bringing the vaccination-reduced reproductive number down to one may not be sufficient to control the disease. Also with some parameter values stability changes can be observed even with forward bifurcation. We would be able to know the minimum portion of susceptible class to be vaccinated in order to have a control over a disease by knowing the bifurcation ii point of backward bifurcation. If so we can avoid a sudden jump in the number of infectives due to an unstable middle branch of bifurcation. i i i Table of Contents Abstract ii Table of Contents iv List of Tables vi List of Figures vii Acknowledgements ix 1 Introduction and Statement of the Problem 1 1.1 General Overview 1 1.2 Research Overview 1 1.3 Thesis Outline 3 1.4 Motivation 4 1.5 Basic Model 5 1.6 Basic SIS model 8 1.7 Basic SIR model 10 1.8 Simplest immigration model without vaccination and the basic reproductive number 14 2 The SVIS model 17 2.1 Introduction 17 2.2 A model of SVIS type without disease-related death 20 2.2.1 Equilibria and stability analysis 23 2.2.2 When do we get a backward bifurcation? 27 2.3 A model of SVIS type with disease-related death, i.e. a ^ 0 29 2.3.1 The case where there are no infective immigrants 32 2.4 Special cases 41 2.4.1 The case without disease fatalities in the absence of immigration, i.e. a = 0 and A = 0 41 iv 2.4.2 The case when all susceptibles are vaccinated immediately and vaccine does not lose its effectiveness i.e. 9 = 0, <f> —> oo 43 2.4.3 The case when the vaccine is utterly ineffective, i.e. a = 1, and there are no disease fatalities, i.e. a = 0 47 2.4.4 The case when the vaccine is absolutely effective a = 0 and no death is caused by the disease, i.e. a = 0 48 3 The SVIR model 51 3.1 The general SVIR Model 51 3.1.1 The vaccine reproductive number, R{<f)) 55 3.2 Endemic equilibria and bifurcations when there are no infective immigrants, p = 0 58 3.3 The general SVIR case: Non-existence of backward bifurcation 62 3.4 Outgoing immigrants and backward bifurcations 63 3.5 Special case: The SVIR model when there are no disease fatalities but with infective immigrants, i.e. p ^ 0 and a = 0 67 3.5.1 Stability analysis 69 4 Conclusion and Discussion 70 A Bifurcation 73 A . l Hopf Bifurcation 73 B Transcritical bifurcation of SVIS model with p = 0 77 C Derivation of Stability Formula of two-dimensional Hopf bifurcation [5] 79 Bibliography 83 v List o f Tables 1.1 E p i d e m i o l o g i c a l of character is t ics of the R N A viruses ( s u m m a r i z e d i n Fields) 6 2.1 S u m m a r y of n o t a t i o n 19 2.2 E x p l i c i t roots of the cub ic p o l y n o m i a l p(X) — A 3 + a \ \ 2 + a 2 \ + 03, w i t h a\, a,2, a n d 03 rea l . a\ = 3a, a2 = 3fr, a = a? — b, f3 = 2 a 3 — 3ab + 03, 4> = ( l / 3 ) s i n - 1 / 3 / [ 2 a ; 3 / 2 ] ( - T T /6 < ^ < TT/6), tp = (1 /3) c o s h " 1 ^ ! / ^ 3 / 2 ] } a n d 6 = (1 /3) s i n h - ^ / t y p t - a ) 3 / 2 ] } , [11] 34 3.1 S u m m a r y of n o t a t i o n 53 v i List of Figures 2.1 The diagram for SVIS model 19 2.2 Bifurcation of SIS model when there is no disease-fatality ,with various frac-tion infective from immigration, p. Prom the top left, various p values are used : p = 0.1,p = 0.2,p = 0.7 and p = 0.3 in clockwise direction respectively. Other parameters are fixed at A = 3, (3 = 1, a = 0.2, p = 0.1, 7 = 12, 9 = 0.5 and A = 2 28 2.3 Bifurcation of SIS model with various fraction infective from immigration, p. From the top left, various p values are used : p — 0.1,p = 0.2,p — 0A,p = 0.48,p = 0.5 and p = 0.3 respectively. Other parameters are fixed at A = 3,-/3 = l , cr = 0.2, p = 0.1, 7 = 12, 6> = 0.5, A = 2 33 2.4 Bifurcation of SIS model when there is no disease-fatality ,with various frac-tion infective from immigration, p. From the top left, various p values are used : a = 0,a = 0.05,a = 0.15 and a — 0.1 in clockwise direction respec-tively. Other parameters are fixed at p = 0, A = 3, (3 = 1, a — 0.2, p = 0.1, 7 = 12, 9 = 0.5 and A = 2 39 2.5 Bifurcation of SIS model in the absence of immigration and when there is no disease fatality, i.e. a = 0, E = 0 42 2.6 Bifurcation when all susceptibles are vaccinated immediately and vaccine does not lose its effectiveness i.e. 9 = 0, <?!>—> 00. Left: the case with no infective immigrant, i.e. p = 0 Right: the case with p = 0.1. Other parameter values are fixed at a = 0.2, p = 0.1, 7 = 12, A = 2, K = 20 while a control parameter is a contact rate, (3 46 vii 2.7 Bifurcation when vaccine is absolutely effective a = 0 and no death is caused by the disease, i.e. a = 0 50 3.1 The diagram for SVIR model 53 3.2 Outgoing immigrants and Hopf bifurcation 1: Other parameters are fixed at A= -5.55, p = 0.50, a = 0.43, (3 = 2.92, a = 0.34, p = 0.68, 7 = 0.53, 6 = 0.36 64 3.3 Outgoing immigrants and Hopf bifurcation 2: Other parameters are fixed at A = -7.65, p = 0.10, a = 0.73, 0 = 3.77, a = 0.18, p = 0.37, 7 = 0.92, 0 = 0.52 65 3.4 Outgoing immigrants and Hopf bifurcation 3: Other parameters are fixed at A = -2.73, p = 0.47, a = 0.90, p* = 6.26, a = 0.06, p = 0.09, 7 = 2.71, 0 = 0.41 66 viii Acknowledgements I would like to thank Fred Brauer, my supervisor, for always being there and providing many suggestions and comments at all stages of this work. I am also thankful to him for his guidance through my early stage of chaos and confusion. Professor Leah Keshet is the one who introduced the subject of mathematical biology to me at the beginning and provided much useful advice. I learned a lot from Professor Yue-Xian Li by taking his class and doing project: I thank him for giving generously of his time to answer my questions. I thank to the department of mathematics and Institute of Applied Mathematics at the University of British Columbia for giving me the opportunity to pursue my graduate studies. Of course, I am grateful to my parents for their patience and love. Without them this work would never have come into existence (literally). Finally, I wish to thank the following: Richard T.(for friendship); Adriana (my office mate who taught many details of software tricks); Tommy (my funny brother who never loses his smile); all my UBC professors; Jeff (for all the good and bad times we had together); Johann Sebastian Bach (and he knows why); and, again, my mother (without you nothing would have been possible). Vancouver, British Columbia Eunha (Alicia) Shim April 20, 2004 ix Chapter 1 Introduction and Statement of the Prob lem 1.1 General Overview In this thesis we demonstrate how a mathematical model can describe epidemiological phe-nomena and how we can use such a model to analyze endemic states and help eradicate diseases. A model was set up using a system of nonlinear ordinary differential equations and this is analyzed mathematically, and simulated using software. Computer software gives a concrete picture of numerical predictions and indicates the direction for mathematical anal-ysis and proof. 1.2 Research Overview Communicable diseases may be introduced into a population through the migration of infective individuals from outside into the host population. In order to eradicate the disease, numerous efforts could be made, such as improved sanitation, antibiotics, and vaccination programs. The model is developed to describe the behavior of a disease when spreading 1 2 through population with immigration and with a vaccination policy in effect. A vaccination program could mean a program of modifying behavior to avoid risk as well as a clinical vaccination. An interesting phenomena in epidemiology is the possibility of multiple steady states and epidemiological mechanisms that produce them. Analysis of the mathematical model shows bifurcation phenomena, usually involving a threshold parameter. The most common parameter is the basic reproductive number, RQ, which is a dimensionless quantity demonstrating an average number of secondary infections caused by an infective individual introduced into a completely susceptible population. In epidemiology Ro is often a threshold parameter which determines whether a disease will invade the host population or not. Such behaviour can be stated as follows in epidemiological terms: if the average number of secondary infections cased by an average infective is less than one, a disease will die out, but if it exceeds one there an epidemic will result[1]. Often this type of phenomenon involves a transcritical bifurcation, which displays an exchange of stability between the disease-free equilibrium and the endemic equilibrium. Without immigration, the disease-free equilibrium exists for all values of Ro while an en-demic equilibrium exists on only one side of the bifurcation point [10]. In the presence of forward bifurcation, an endemic equilibrium exists only when Ro > 1 in the absence of infective immigration. On the other hand, when backward bifurcation is present, an endemic equilibrium also exists when Ro < 1, so there can be an endemic states even when Ro < 1 for some parameter sets. In a system with a backward bifurcation, when Ro < 1 but is still close enough to one, there exist both stable and unstable endemic 3 equilibria, while for Ro far below one there exists only stable disease-free equilibria. Thus, increasing Ro to just above one will create a dramatic jump in the number of infectives [10]. In this case, in order to eradicate the disease, Ro needs to be reduced far below than one until it passes a so-called saddle-node bifurcation at Ro < 1, [10]. The main aims of my thesis are: • To provide examples of backward bifurcation of an epidemic model and to explain its effect on vaccination policy; • To show how each different model of different disease types (mainly SIS and SIR types are considered) can produce significantly different bifurcation types; and • To demonstrate how standard and bilinear incidence can produce different types of numerical results. In achieving these aims, I use three different models and show numerical and analytic result from each case. The computer packages I used are M A T L A B (for random search of appropriate parameter sets), XPPAUT, and A U T O . 1.3 Thesis Outline Chapter 1 covers preliminary examples that will serve as an introduction to some termi-nology and results already well known. I begin Chapter 2 by presenting an epidemiological model of the SVIS type when the immigration factor is included. A disease is classified as an SIS type when it confers no immunity - one can get infected later although he (or she) 4 has recovered from the infection previously. In this model the total population is divided into three groups - susceptibles, infectives, and vaccinated. Also there is a constant flow of immigrants from the outside into the host population. Chapter 3 covers the epidemiological model of the SVIR type when the immigration factor is considered. A disease is classified as an SIR type when it confers immunity - one will not be infected again once he (or she) recovers from infection. We focus on the different behaviour of the bifurcation diagram from the SVIS case. Various numerical results are provided. 1.4 Motivation In 1665-1666 London suffered two epidemic disasters. In the spring and summer of 1665 an outbreak of the plague that is now identified as Bubonic Plague spread out until thousands had died. The plague indeed started in the East and invaded through Europe quickly. Bubonic Plague was called the Black Death since the victim's skin became black in patches when infected. The symptom include inflamed glands, vomiting, swollen tongue and split-ting headaches. The spread of disease started slowly at first but 43 had died by May 1665. Then 6137 people died in June, 17036 in July and at its peak in August, 31159 people died. In all the total population was reduced by 15%. The characteristic features of the Great Plague were its short incubation period and the sudden emergence of the epidemic. In fact incubation took less than a week. Since the Great Plague many scientists have paid attention to models for epidemics in mathematical terms. Also the epidemic cholera was originally isolated to India mainly, and outbreaks occurred recurrently in cyclic fashion. However, through the centuries, the range 5 of cholera infection spread rapidly. By the 18th and 19th centuries, in Moscow (1830), and North America (1832), cholera epidemics were reported. Some epidemic diseases such as measles and influenza can be prevented by vaccination and many researchers have tried to study the periodic behavior of measles for instance. Table 1.1 summarizes characterizes of some disease whose infection can be prevented by vaccination. 1.5 Basic Model To begin this section I will present one of most simplest epidemiological models of SIR type. This model is based on a few assumptions below: • An average infective individual make an appropriate contact sufficient to transmit infection with /3N others per unit time. • The probability that a contact of an infective individual is with a susceptible is S/N. • A fraction 7 of infectives recover per unit time. • The demographic effect can be ignored since the epidemic time scale is very short. Based on these assumptions Kermack and McKendrick proposed the following model in 1927: S' = -psi I' = psi - il Rl = 7 / 6 Virus Measles Influenza A Influenza B Influenza C (human) (equine) (avian) Main Species Human Human Equines All avian Infected species Transmission Direct Direct Direct Direct respiratory respiratory respiratory respiratory and oral Incubation 8-12 days 2-3 days 1-2 days -Period Infectious 5-8 days < 14dd in 7dd in large Months (in Period (in young amounts, wild ducks) naive hosts) < 14dd in small amounts Infectious N / A l-3dd 1-5 days N / A Period Extent of Lifelong after Prolonged, but Clinical Highly variable immunity natural strain immunity for (rapid clearance following infection dependent >12 mo, but in some infection/ replication and species, vaccination transmission occurs < 12 mo prolonged enteric infection Current Infant Vaccination Vaccination Vaccination control Vaccination of elderly quarantine quarantine strategy and at risk Min 300,000 Probably much As for human "Small" population size for larger than measles for a (persistent infectious) persistence given strain Epidemic Seasonal and Seasonal Seasonal -dynamics multi-annual epidemics epidemics with waves epidemics Table 1.1: Epidemiological of characteristics of the RNA viruses (summarized in Fields) 7 In this model R can be determined once / is determined since the differential equation for R only depends on I. Hence we consider the equation for S and I first as below: S' = -0SI I' = (0S - 7 ) / Since solving these equations analytically would be rather difficult one can use a qualitative approach to study the behavior of the system. This approach is well presented in [1]. We begin by dividing the equation of I' by the one of S': T _dl _ (pS - a)I _ a_ S1 ~ dS ~ -pSI ~ + PS By integrating both sides we obtain I = -S+^\ogS + c and we define the function V(S,I) as where Now it follows that V(S,I) = S + I-^logS a K — S0 OL CX CX CX I m a x = S0 + h - log So ~ -jj + ~jj log -Some diseases such as measles confer immunity so recovery from such diseases provides immunity against reinfection. We shall examine a general SIR model with births and deaths 8 for a disease S' = fiK- BSI - pS I' = pSI-(ji + a + y)I (1.5.1) R' = 7 J - ^i? with 5 + I + R = K. \x + 7 + a boo 1.6 B a s i c S I S m o d e l Some diseases confer immunity so recovered individuals gain immunity against the disease. This type of disease can be modelled by SIS type where S and I stand for susceptibles and infectives respectively. Here we proposed a simple model of SIS type that includes births in the susceptible class and deaths from two classes with the rate proportional to each class. Our model is S' = -pSI + n(K -S) + jI r = PSI-(a + fi + 7 ) / where the total population size, N, is defined as sum of population in two classes, i.e.N(t) — S(t)+I(t). We assume that a disease may be fatal to some infectives so deaths due to disease can be included in a model using the disease-related death rate from infective class, a. The recovered individuals are removed from the infective class at the rate, 7, and they don't acquire immunity. The birth rate is constant and proportional to the carrying capacity, K, and the natural death rate is fi. We can start a qualitative approach by studying steady states. The equilibrium condition can be found by setting the right hand side of equation 9 to be zero: -BSI + u.{K-S) + jI = 0 pSI-(a + fi + 7)/ = 0 The second equilibrium condition factors, and we get one disease-free equilibrium ,1 = 0 and the other endemic one that satisfies PS — (a + ^  + 7) =0 assuming (3S — (a + ^  + 7) < (3K. Hence by solving the two equilibrium conditions one can find the two following steady states, a disease-free and an endemic one: (I*,ST) = (0,K) P{a + u- + 1) ' } To study the local stability of a fixed point we linearize our system first. The Jacobian matrix is given by j = ( - u - 3 I -BS \ fil B'S-i-y + fi + a) The Jacobian matrix at the disease-free equilibrium is given by -H -PK + 7 0 pK - (7 + a + n) By looking at eigenvalues, one can easily see that the disease free equilibrium is stable if PK < (7 + fj, + a). Now we can turn to an endemic equilibrium and study about its stability. The Jacobian matrix evaluated at an endemic equilibrium is given by - r ^ y - ( " + « ) We can easily see that the trace of this matrix is negative and the determinant is positive as long as PK — (7 + ^  + a) > 0. 10 • URo < 1 — the disease-free equilibrium of 1.7.1 is stable, and — the endemic equilibrium does not exists. • If #0 > 1 — the disease-free equilibrium of 1.7.1 is unstable, and — the endemic equilibrium exists, and is asymptotically stable. where the basic reproductive number is defined as Ro = p+^+a- Note that the basic reproductive number depends on the contact rate /3, carrying capacity, K, and the rates at which an infective is removed from its class. 1.7 B a s i c S I R m o d e l Some diseases confers immunity so recovered individuals gain immunity against the disease. This type of disease can be modelled by SIR type where 5, I, and R stands for susceptibles, infectives and the recovered respectively. The first SIR type model was proposed by H.E. Soper. He assumed that the total population size is constant and that the birth rate and death rate are constant. His model is S' = -0SI + \±K r = psi - 7/ R' = 7 / - pK where P is contact rate and 7 is recovery rate. However this model was criticized since there is no proper reason why one should link deaths from recovered class to births from 11 susceptibles, [1]. Also R(t) can become negative when R(0) and 1(0) are sufficiently small, which is another weakness of Soper's model. In 1932, Kermack and McKendrick proposed another model of SIR type that includes births in the susceptible class and deaths from all three class with the rate proportional to each class. Their model is S' = -BSI + n(K ~ S) I' = BSI -yl-iil R' = 7 J - nR where the total population size, N, is defined as the sum of the populations in the three classes, i.e.N(t) — S(t) + I(t) + R(t). Their model treats only the case where there is no death due to disease. However one can assume that a disease may be fatal to some infectives so deaths due to disease can be included in the / equation. Thus in addition to the demographic effects we can also consider removal of some infectives via disease-related death. Based on this assumption we can consider a more general model that was proposed by Hethcote in 1976, [6]: S' = \±K- BSI - fiS I1 = BSI — ( 7 + /z + a)I (1.7.1) Rl = 7 J - nR Here we assume that a is the fraction of infectives dying from infection and 7 is the rate at which infectives recover with acquired immunity. The natural death rate is fi and the birth rate, fiK, is assumed to be constant. Let us start to analyze the model qualitatively 12 by eliminating the R variable from the system and reducing to the two-dimensional one: S' = -0SI + p{K - S) I1 = 0SI-[n + ii + a)I We can ignore the equation for R since R can be determined once S and I are known [1]. We can start the qualitative approach by studying steady states. The equilibrium condition can be found by setting the right hand side of the equations to be zero: -0SI + n{K - S) = 0 0SI - (7 + » + a)I = 0 The second equilibrium condition factors, and we get one disease-free equilibrium, 1 = 0, and the other endemic one that satisfies (3S— ('j+p+a) = 0 assuming (3S = (y+p+a) < 0K. Hence by solving the two equilibrium conditions one can find the following two steady states, a disease-free and an endemic one: (I*,S*) = (0,K) ,T* c*\ .. / M / ^ ~ (7 + M + a)] 7 + M + /3(7 + M + «) ' 0 ' To study the local stability of a fixed point we linearize our system first. The Jacobian matrix is given by J = ( -U.-0I -0S \ 01 / 3 5 - ( 7 + M + «) The Jacobian matrix at a disease-free equilibrium is given by -IM -0K 0 0K - (7 + /J, + a) 13 By looking at eigenvalues, one can easily see that the disease free equilibrium is stable if QK < (7 + n + a). Now we can turn to an endemic equilibrium and study its stability closely. The Jacobian matrix evaluated at an endemic equilibrium is given by The trace of this matrix is always negative and the determinant is positive. as long as BK — (7 + fi + a) > 0, that is, the same condition as the one for existence of an endemic equilibrium. In summary there is threshold behavior in this model: • liRo < 1 — the disease-free equilibrium of 1.7.1 is stable, and — the endemic equilibrium does not exists. • If iJo > 1 — the disease-free equilibrium of 1.7.1 is unstable, and — the endemic equilibrium exists, and is asymptotically stable. where the basic reproductive number is defined as RQ — _ ^ _ a - Note that the basic reproductive number depends on the contact rate, 0, carrying capacity, K, and the rates at which an infective is removed from its class. 14 1.8 S i m p l e s t i m m i g r a t i o n m o d e l w i t h o u t v a c c i n a t i o n a n d t h e b a s i c r e p r o d u c t i v e n u m b e r Let us now consider an SIS type disease in the absence of vaccination but with a constant flow of incoming immigrants. Let S(t) be the number of members of the population who are susceptible to an infection at time t and let I(t) be the number of members who are infective at time t. The total population size at time t is denoted by N(t), with N = S + I for the SIS model(the disease confers no immunity). When an infective make contact, the probability of producing a new infection is since a new infection can be made only when a contact is made with a susceptible. Thus the rate of producing new infections is J3N • § • I = /3SI, where (3 is a constant. For simplicity, we assume that there is no disease related death. The population is replenished via two ways, the birth and immigration. We assume that newborns enter the susceptible class at a constant rate and there is a constant incoming flow of immigrants where some portion, p, of immigrants are infectives. In summary, the assumption we have in this model is as following: • S, I and N are the number of susceptibles, infectives and total population, respec-tively. • There is a constant flow of A new members into a population in unit time, where a fraction p > 0 of A is infective (0 < p < 1). • There is a constant per capita natural death rate p > 0 in each class. • A fraction 7 > 0 of infectives recovers in unit time. 15 • 0 is the infectious contact rate per person in unit time • A is the constant natural birth rate; all newborns come into the susceptible class We can now formulate this model, dividing the total population into two classes - suscep-tibles (S) and infectives(I): S1 = (1 - p) A + A - BSI - LiS + 7 / (1.8.1) I' = pA + BSI + 7 ) / Note that the total population is the sum of two classes, susceptibles and infectives, thus N{t) = S(t) + I{t) It follows that N' = S' + I' = A + A - fiN Since lim^oo N(t) = K, the model is asymptotically stable and we can reduce the size of model by replacing S by K — I, [4]. The model is now equivalent to I'=pA + fi{K - 1)1 - (jx + 7 ) J (1.8.2) = -BI2 + (BK-u--j)I+pA By setting I' = 0, we can solve for the two equilibria when {0K — fx — 7 ) 2 — 4(—B){pA) > 0: _ (BK - Li - 7) ± y/{0K - n - 7 ) 2 - 4(-/?)(M) h ' 2 ~ 2B where I± > 0 and i | < 0. Let us consider an endemic equilibrium, I±. We can rewrite II by letting * = BK - fi - 7, [3]: 16 lim I\ * + |*| = < 0, if * < 0; f , i f * > 0 (1.8.3) P ^ O * 2/3 From this expression we obtain a threshold parameter, and deduce the basic repro-ductive number, RQ: Ro = Using RQ, we can rewrite the equilibria. 0K „ (Ro - + 7) ± V(Ro - 1)2(M + if + 4/3pA Note that > 0 and J | < 0 for all i? 0 > 0. Chapter 2 The SVIS model 2.1 In t roduc t i on Let us now consider an SIS type disease when a vaccination program is in effect and there is a constant flow of incoming immigrants. Let S(t) be the number of population members who are susceptible to an infection at time t, I(t) be the number of members who are infective at time t, and V(t) be the number of members who are vaccinated at time t. The total population size at time t is denoted by N(t), with N = S + 1 + V for the SIS model (the disease confers no immunity). Assume that each infective makes /3N contacts sufficient to transmit infection in unit time, where (3 is a const. When an infective makes contact, the probability of producing a new infection is -j|, since the new infection can be made only when a contact is made with a susceptible. Thus, the rate of producing new infections is (3N • jj • I = (3SI. The susceptible population is vaccinated at a constant rate, (f>, and the rate at which the vaccine wears off is 9. We assume that there can be disease related deaths as well as natural deaths unrelated to the disease. The population is replenished in two ways, birth and immigration. We assume that all newborns enter the susceptible class at a 17 18 constant rate, A, and that there is a constant incoming flow (A) of immigrants where some portion of immigrants, p, is infective. In summary, the assumptions we have in this model is as follows: • S(t), I(t), V(t), and N(t) are the numbers of susceptibles, infectives, vaccinated, and total population at time £, respectively. • There is a constant flow of A new members into the population per unit of time, where fraction p of immigrants is infective (0 < p < 1). • The vaccine has the effect of reducing infection by a factor of cr, so that a = 0 means that the vaccine is completely effective in preventing infection, while a = 1 means that the vaccine is utterly ineffective. • The rate at which the susceptible population is vaccinated is <f>, and the rate at which the vaccine wears off is 6. • The disease can be fatal to some infectives and we define a to be the rate of disease-related death. • There is a constant per capita natural death rate p > 0 in each class. • Fraction 7 > 0 of infectives recovers in unit time. • PN is the infectious contact rate per person in unit time. • A is the constant natural birth rate, with all newborns coming into the susceptible class. S(t) Number of susceptibles at time t V(t) Number of vaccinated individuals at time t I{t) Number of infectives at time t R(t) Number of recovered people with immunity at time (t) A Number of immigrants P Portion of infectives among immigrants A Birth rate 0 Contact rate 7 • Recovery rate <t> Vaccination rate a Factor by which the vaccine reduces infection e Rate at which the vaccine wears off Natural death rate unrelated to the disease a Disease related death rate Ro Basic reproductive number R<t> Vaccination reproductive number Table 2.1: Summary of notation SSI pA S MS i oS 1 5~ ( a or V i '"T'T fit al Figure 2.1: The diagram for SVIS model 20 We can now formulate this model, by dividing the population into three classes - sus-ceptibles (S), vaccinated (V), and infectives (I): S' = (1 - p)A + A - PSI - (/x + 4,)S + ~fI + 6V I'=pA + PSI + apVI -{p + j + a)I V' = <f>S- apVI -(p + 0)V (2.1.1) Note that the total population is the sum of three classes: susceptibles, infectives, and vaccinated: N(t) = S(t) + I(t) + V(t) Thus, it follows that N' = S' + I' + V = A + A-pN-al The equation 2.1.1 is the SVIS model we will use to investigate the behavior of an SIS type disease throughout this chapter. We will study first the simpler case of no disease-related death. 2.2 A m o d e l of S V I S type wi thou t disease-related dea th If the disease causes no fatality (a = 0), then we can simplify the model 2.1.1 by letting a = 0 as follows: S' = ( i - p)A + A -psi -{p + cf>)S + 1 i + ev I' = pA -f PSI + GPVI - (p + 7 )7 (2.2.1 V = (f>S -apVI-{p + 6)V 21 with N' = S' + I' + V = A + A - LIN The system is now asymptotically autonomous, so we can define the limit value of iV as follows: lim N(t) = = K Definition 2 .2 .1 . Consider the differential equation § = * = fM (2-2.2) where x = x(t) G M™ is a vector valued function of an independent variable (usually time) and /:[/—> IRn is a smooth function defined on some subset U CW1. Systems of the form 2.2.2, in which the vector field does not contain time explicitly, are called autonomous, [12]. Definition 2.2.2. The system y' = f(t,y) is said to be asymptotically autonomous on the set if and only if 1. limt_n3o/(t,y) = h(y) for y E fi and this convergence is uniform for y in closed bounded subsets of fl. 2. For every e > 0 and every j / e f i there exists a 6(e, y) > 0 such that |/(i , x) — f(t, y)\ < e, whenever \x — y\ < 5 for 0 < t < oo, [2]. 22 According to the theory of asymptotically autonomous system we can reduce it to two-dimensional system by letting replacing S with K — I — V: I' = PA + 0[K - I - (1 - a)V]I -{p + 7) J V' = 4>[K - I - V] - aPVI -{p + 9)V (2.2.3) This is the system that we will analyze in order to find the basic reproductive number. At an equilibrium, we can set the right hand side of equations 2.2.3 to be zero, which gives the equilibrium conditions: P(l - cr)V = ?j + P(K-I)~(p + 7) 4>{K -I) = V\aPI +{p + 0+ (/))} One can solve for V from the second condition and substitute V into the first one: pA p(l - o-)4>{K -I) = [aPI +(p + 6 + cf>)}[^r+ P{K -I)-{p + 7)] (2.2.4) P<j>(l -a){K- 1)1 = [api + {p + 0 + <f>)]\pA + pI{K - I ) - I{p + 7)] After complicated computations we omitted here, the endemic equilibrium condition, 2.2.4, becomes a cubic polynomial of / as below: /(/) = EIZ + BI2 + CI + D = 0 (2.2.5) with E =pa B ={p + 0 + a<t>) + a(p + 7) - <rpK C = -paA-K{p + 0 + a<l)) + (p + j){p + 9 + <f>)/p pA(p + 0 + 4>) D = -P 23 P r o p o s i t i o n 2.2.1. If the system 2.2.1 has three distinct endemic steady states then two conditions on the parameters must be satisfied: ((j. + 6 + a(j>) + a{n + 7) - crBK < 0 and - paA - K(fi + 0 + a<j>) + (// + + 6 + <j))/B > 0. Proof. Consider the equilibrium condition 2.2.5. Note that the sign of some coefficients of the polynomial are fixed for all nonnegative parameters: E > 0 and D < 0 for all nonnegative parameters. If D ^ 0, then there are either one or three positive roots since /(0) < 0 and lim^oo / ( / ) = oo. Through differentiating 2.2.5 with respect to / , we get /'(/) = 3EI2 + 2BI + C. If /(/) = 0 has three positive roots, then it follows that /'(/) = 0 has two positive roots according to Rolle's theorem. Thus, B < 0 and C > 0 are the necessary conditions for the system 2.2.1 to have three endemic steady states. • 2.2.1 Equilibria and stability analysis P r o p o s i t i o n 2.2.2. The sign of the slope of the bifurcation curve, I vs P, of the system 2.2.1 is determined by the sign of (3EI2 + 2BI + C) where . E=pa B=(n + 9 + a4>) + O{LI + 7) - o-BK C = -paA-K(Li + 6 + o-<t>) + (Li + 7){Li + e + (f>)/p The slope of the bifurcation curve is positive when (3EI2 + 2BI + C) > 0 and negative when [3EI2 + 2BI + C) < 0. 24 Proof. Recall that the equilibrium condition of system 2.2.3 is the cubic equation /(/), 2.2.5. We have already shown that B < 0 and C > 0 is a necessary condition that /(/) has three positive equilibria using Rolle's theorem. Now consider a bifurcation diagram using (3 as a control parameter. Implicit differentiation of 2.2.5 with respect to f3 gives: / „ ^ r 2 ~r,T ^ d l TodE r2dB TdC dD (3EI2 + 2BI + C ) T r - 1 ^ - I 2 ^ - 1 ^ - ^ = - / V + I\*K) + I f (^ + 7)(M + g + 0) \ + pAfr + e + fi (32 J /3 2 = a2I(K-I) + ^ ± ^ [ I ( p + 1)+pA} (2.2.6) We can easily see that the right hand side of implicit differentiation 2.2.6 is positive, since K > I. Thus, the sign of slope of the bifurcation curve, is determined by the sign of (3EI2 + 2BI + C). • The slope of the bifurcation curve of I using (3 as a control parameter is infinitely large when dl From 2.2.6 the slope becomes infinitely large when 3 £ I 2 + 2BI + C = 0. (2.2.7) where E = /3a B = {p + 0 + cr</>) + (T(JM + 7) - o~(3K C = -pa A -K(p + 0 + a<f>) + {p + 7 ) (M + 0 + <fi)/(3 25 It is easily seen that E > 0 for positive parameters. Hence one can easily see that there exist two distinct positive roots of equation 2 . 2 . 7 i f i ? < 0 , C > 0 and y/B2 — 3EC > 0: I„ = -_B±VWEMc (2 2 8 ) These two roots are indeed bifurcation points of 2.2.3 at which the slope of bifurcation becomes infinitely large. One can also study the stability of equilibrium from the linearization of 2.2.3: J = ( -2[3I- (1 -a)pV-(u. + "f) + 0K -(l-a)pi \ y -(<f> + a(3V) -(Li + 8 + cf> + o-l3I) ) Using the equilibrium condition, = 0K — 01 — 0(1 — a)V — (LI + 7 ) , one can evaluate the Jacobian matrix at equilibrium: J \ , - ( ) P.2.10) y-(cf> + a0V) -(Li + 9 + cf> + a0I) J Note that this matrix has negative trace. The determinant of this matrix is det(J) = a(0I)2 + 3I(n + 9 + <j>) + ^ (LL + 9 + (j) + afil) - (1 - o)<j>BI - (1 - a)cj02VI = QI\oQI +(fi + 9 + acf>)- 0a(K -I) + (LI + i)<r\ + ^ (fi + 8 + <j>) Let us further simplify the determinant using equilibrium condition, ^-+B(K- I)a - BVa(l -a)-(u. + j)a = 0. 26 det = BI[cBI + {p + 6 + <r<f>) - Bo(K -I) + {p + 7)0-] + ^ - Q u + 9 + <f>) = B\2oBI + (p + 9 + acj)) + a{p + 7) - <TBK\ + ^ ( M + B + 4>) = (3I[2EI + B] + ^ ( t i + 9 + (f>) (2.2.11) = 8I[2EI + B - jf] = j[2EIs + BI2 - D] = ^[I2(2EI + B)-D] Recall that D < 0 and E > 0 for all positive parameters. Thus the determinant is positive if B > 0. Since the trace is negative, the steady states are asymptotically stable when B > 0 where B = (fj, + 9 + a<j)) + a(p + 7) - aBK. Proposition 2.2.3. The exchange of stability occurs when the slope of bifurcation curve, 8 vs I, changes. Proof. Let us recall the equilibrium condition, 2.2.5. Throughout a stability analysis from above, we know that the determinant of zero indicates the threshold for stability changes: Second, when the sign of the slope of the bifurcation curve changes, one can find a threshold point by letting the left hand side of equation 2.2.6 be zero and solving for I: 2EIi + BI2-D = 0 (2.2.12) 3EI2 + 2BI + C = 0 3EI3 + 2BI2 + CI = 0 (2.2.13) 27 Now, we can show that these two conditions are equivalent by using the equilibrium condi-tion EI3 + BI2 + CI + D = 0. {3EI3 + 2BI2 + CI) - (2EI3 + BI2 - D) = EI3 + BI2 + CI + D = 0 (2.2.14) which is the equilibrium condition. In summary, using the equilibrium condition, EI3 + BI2 + CI + D = 0, one can show that when the sign of the slope of the bifurcation curve changes, its stability changes too. • 2.2.2 When do we get a backward bifurcation? The idea underlying this section is that an infinite slope of the bifurcation curve signals a backward bifurcation. The slope of the bifurcation curve is infinite when /(/) = 3EI2 + 2BI + C = 0. Thus if /(/) has two positive distinct real roots, we get backward bifurcation. The sufficient conditions for /(/) to have two positive real roots is as follows; /(0) > 0, ^ > 0 and B2 - 3EC > 0 3E which is equivalent to P< < ti{fi + j){p + 0 + <f>)-B(A + A ) ( l + <r<t>) ap(3A a/3(A + A ) - ap(p + 7) - p(p + 0) ap and p > 1 3a2p/3A 3(3a(A + A) {p + 0 + a(f>) + 3pa(p + j){p + 9 + a<p) - p[(p + 0 + a<l>)+a{p + 'y) af3{A + A) 2l (2.2.15) Prom Figure 2.2, one can easily notice that we lose backward bifurcation as we increase p when other parameters are fixed. 28 Figure 2.2: Bifurcation of SIS model when there is no disease-fatality ,with various fraction -infective from immigration, p. From the top left, various p values are used : p = 0.1,p = 0.2,p = 0.7 and p = 0.3 in clockwise direction respectively. Other parameters are fixed at A = 3, 0 = 1, a = 0.2, \i = 0.1, 7 = 12, 9 = 0.5 and A— 2. 29 2.3 A m o d e l o f S V I S t y p e w i t h d i s e a s e - r e l a t e d d e a t h , i .e. a ^ 0. We now consider the original model which includes non-zero disease-fatality. Recalling the original system of differential equation, let us remind the system of differential equations, 2.1.1: S' = (1 -p)A + A - BSI - (/i + 4>)S + + I' = pA + BSI + aBVI - (^  + 7 + a)I V = <l>S-a8VI-(fi + e)V Recall that the total population is the sum of three classes, susceptibles, infectives and vaccinated, i.e. N(t) = S(t) + I(t) + V(t) Thus it follows that N' = S' + I' + V = A + A-LiN-al We can get an alternate but yet equivalent model by replacing S by N — V — I. Now the model is: I' = pA + BI[N - I - { I - a)V] -{Li + j + a)I V' = c/)(N - I ) - aBVI -(ii + 0 + (f))V N' = A + A — LIN — al (2.3.1) We can write the equilibrium conditions by letting the right-hand side of 2.3.1 be zero. The 30 equilibrium conditions are pA + 0I[N - I - (1 - a)V] ={LI + 7 + a)I <f>(N-I) ={n + 6 + 4> + a0I)V LIN + al =A + A From the second and the third equilibrium conditions we can solve for N and V: A + A-al N V = <f>(N -1) H + 9 + 4> + a 01 <f>[A + A - (a + n)I\ tl(Ll + 6 + (f> + (701) Eliminating N and V by substitution of these expressions into the first condition and simplifying, we obtain an expression involving I for the equilibrium condition of the form pA + 01 A + A-al (j>[A + A - (a + LI)I) [ a ) Li^ + e + ^ + aPI) = (n + 7 + a)I We can further simplify by multiplying + 9 + <f) + (701) in order to get a cubic equation for the equilibrium values of / of the form EI3 + BI2 + CI + D = 0 where E = 0(7 [a + n) B = -(A + A)cr/3 + (a + LL){LI + 9 + a(f>) + + 7 + a) + 7 + a) {fj, + 9 + (j>) C = -LipaA - (A + A)(/j + 8 + a4>) + -pAn{n + 9 + (f>) 0 D = 31 In order to study the stability of steady states we start a qualitative approach by lineariza-tion of 2.3.1: / B[N -(l-a)V-I]-pi-{fx + j + a) -B(\ - a)I [31 \ J = -(4> + a BV) V -a -o-f3I -(n + 6 + 4>) 4> 0 -p J Using the equilibrium condition, pA + BI[N - I —{I — a)V] = (p + j + a)I, we can rewrite the matrix as: ( -B(l-a)I (31 \ J\i* = -(</> +<rpV) -aj3I-{p + e + <l>) 4> \ -a 0 -p j After a complicated computation, we can obtain its characteristic equation involving I in it: where X3 + ai\2 + a2X + a3 = 0 ai = (1 + a)BI + ^ + 2p + 6 + (f> a2 = (BI + ^)(p, + e + (f> + (T0I) + (p. + e + (l> + aBI)p + p(J3I + ^) - BI(1 - a)(<f> + aBV) + aBI a 3 = n(J3I + ^)(p + 9 + cf> + aPI) - p(3I(l - a)(cf> + a0V) + a[3I{p + 6 + a(/) + a (31) By the Routh-Hurwitz Criterion, the steady state is globally stable if and only if a\ > 0 a$ > 0 and a\a2 > 03. The Figure 2.3 demonstrates a few cases where an equilibrium graph loses its stability as the vaccination rate, </>, increases and becomes stable again. At the point where it loses local 32 stability first, Hopf-bifurcation occurs and a periodic solution appears for some <f> values. Table 2.2 summarizes the explicit roots of the cubic polynomial, A 3 + AX2 + BX + C — 0 with real a\, a2 and a^. Also the theorem about the normal form of Hopf-bifurcation and the method to detect various bifurcations in 1-dimensional space is included in Appendix. Theorem 2.3.1. Routh-Hurwitz Theorem, [2] Consider the characteristic equation \XI - A\ = Xn + M " - 1 + • • • + ifc-iA + bn = 0 determining the n eigenvalues X of n x n square matrix A, where I is the identity matrix. Then the eigenvalues X all have negative real parts if where A i > 0, A 2 > o, • - , A „ > 0 1 0 0 0 0 . . 0 bz b2 h 1 0 0 . . 0 h b4 h b2 h 1 . . 0 y b2k-i b2k_2 b2ks b2k-4 b2k-5 ^2fc-6 • • • bk J 2.3.1 The case where there are no infective immigrants It is worthwhile to consider the case without infective immigrants since in this case the system will have a disease-free steady state that would not exist otherwise. This model was proposed by Kribs-Zaleta and Vekasco-Hernandezin, [10]. If there is no infective portion 33 Figure 2.3: Bifurcation of SIS model with various fraction infective from immigration, p. From the top left, various p values are used : p — 0.1,p = 0.2,p — 0.4,p = 0.48,p = 0.5 and p = 0.3 respectively. Other parameters are fixed at A = 3, 0 — 1, a = 0.2, ii = 0.1, 7 = 12, 0 = 0.5, A = 2 34 a > 0 0 = 0 Ai A 2 A 3 = —a = (3a) 1/ 2 - a = - ( 3 a ) 1 / 2 - a < 2a 3 / 2 Ai A 2 A 3 = 2a 1 / 2 sm<f> — a = - 2 a 1 / 2 sin(7r/3 + 4>) - a = 2a1>2 sin(7r/3 - <f>) - a /3 > 2a 3 / 2 Ai A 2 A 3 = —2a1/2 coshi/) — a = a 1 / 2 coshi/' — a + z^a) 1 / 2 sinh V = a 1 / 2 cosh^ > — a — i ^ a ) 1 / 2 sinhi/> 0 < - 2 a 3 / 2 Ai A 2 A 3 = 2a 1 / 2 cosh'i/; — a = —a 1 / 2 cosh'*/' — a + i(3a) 1/ 2 siring = —a 1 / 2 cosh^ — a — i(3a) J/ 2 sinh^ a = 0 —oo < 0 < oo Ai A 2 A 3 = - 0 1 / 3 - a = / 3 1 / 3 / 2 - a + 3i/32/3/4 = 01/3/2-a-3i02'3/4 a < 0 —oo < 0 < oo Ai A 2 A 3 = - 2 ( - a ) 1 / 2 s i n h ^ - a = ( -a) 1 / 2 s inh6'-o + i(-3a)1/2cosh6> = ( -a) 1 / 2 sinh 9 - a - i{-3a)V2 cosh (9 Table 2.2: Explicit roots of the cubic polynomial p(X) = A 3 + aiA 2 + a 2A + a 3 , with ai, a 2 , and a 3 real, ai = 3a, a 2 = '36, a = a 2 - 6, /3 = 2a3 - 3a6 + a 3 , 0 = (1/3) s i n - 1 /3/[2a3/2] ( - T T / 6 < <f> < TT/6), if> = (1/3) cosh"1 {|/3|/[2a3/2]} and 9 = (1/3) sinh"1 {/3/[2(-a)3/2]}, [11] 35 from immigrants, i.e. p = 0, then our equation becomes: S' = A + A - asi -{p + <j>)S + + I' = [3SI + oBVI - (p + j + a) I V = <j)S - oBVI -(p + 9)V The total population is defined as sum of three classes: susceptibles, infectives and vacci-nated, thus, N(t) = S(t) + I(t) + V(t) Thus it follows that N'(t) = S'(t) + I'{t) + V'(t) = A + A - pN(t) - al(t) As before we can make a similar transformation by replacing S(t) by N(t) — V(t) — I(t). Now the model is: = BI[N --I-(l-a)V}- (p + 'y + a)I V = <t>{N- I) - oQVI - {p f 9 + <j))V N' = A + A - pN - al (2.3.2) We can write the equilibrium conditions by letting the right-hand side of 2.3.2 to be zero. The equilibrium conditions are BI[N - I - (1 - <r)V] ={p + 7 + a)I 4>{N-I) ={p + 9 + cf> + apI)V pN + al =A + A 36 From the second and the third equilibrium conditions we can solve for N(t) and V(t): A + A - al(t) N(t) = V(t) J>(N(t) - I(t)) LL + 9 + 4> + a0i(t) <f>[A + A - ( a + u)I] Ll{Ll + 0 + <f> + (T0I{t)) Eliminating N(t) and V(t) by substitution of these expressions into the first condition and simplifying, we obtain an expression involving I(t) for the equilibrium condition of the form 01 A + A - a l T , <j>[A + A - ( a + u)I] 1 — i — cr = (/i + 7 + a)I Ll Lt(Ll + 6 + <f> + a0I) We can further simplify by multiplying + 9 + <j) + oBT) and factoring out a disease-free equilibrium I* = 0 in order to obtain an endemic condition as the quadratic equation for the equilibrium values of I of the form EI2 + BI + C = 0 where E = fio(a + LI) B = -(A + A)o~B + ( « + + 9 + a(f)) + GLI(LI + 7 + « ) LI{LI + -y + a)(Li + 9 + (p) C = -(A + A)(n + 9 + a(j)) + 0 In order to study the stability of steady states we linearize 2.3.2, obtaining the Jacobian matrix: / 0[N - (1 - a)V - /] - 01 - Gu + 7 + ex) J = -(</> + a0V) —a -0(1 - a)I BI \ -a0I -(Li + 9 + <t>) 4> 0 -Ll J 37 By evaluating at the disease-free equilibrium, I* = 0, we can obtain three real eigenvalues: 8(ii + 9 + o-(f>)(A + A) A l /i(ii + 9 + 4>) A 2 =-{<t> + e +u) - (7 + a + p) A 3 = -\i Note that the first eigenvalue is negative when R(4>) < 1 where vaccination reproductive number is defined by R(<j>) f3{n + e + a<t>){A + A) li(li + 9 + <f))(j + a + /i) Thus the disease-free equilibrium is asymptotically stable if and only of R((f>) < 1. Also by using the endemic equilibrium condition, B[N — I— (1 — u)V\ — (/j. + 'y + a), we can evaluate the Jacobian matrix at endemic equilibria: -dl -3(l-o-)I 01 \ -(<!> +aBV) -o3I-{ii + Q + <t>) 4> -a 0 — fi J with a characteristic equation: A 3 + a iA 2 - r -a2A - t - a 3 =0 where ai = (1 + a) 31 + 2p + 9 + (f> a 2 = (31 + u)(ti + 6 + <l> + o-dl) + 3(ti + a)I - 8(1 - a)(<f> + aSV) a 3 = ii3I\n + 0 + o-cf) + o-pi + aP(a - 1)V] + a3I(ji + 9 + acf) + a81) For this SVIS model there is a transcritical bifurcation at (u. + 9){u.(i + a + u)-B(A + A)} a8(A + A) — /n(7 + a + ji) 38 and this is demonstrated in Figure 2.4 below. One can easily see that the lower branch of the bifurcation curve is negative for <f> < ^+fplA+lyi4^a+X"^' a n c * c ° i n c i d e s w i * n t n e disease-free equilibrium value of zero at <f> = ^+fp\A+l^-^^a+^)A^ • Also t n e disease-free equilib-rium is locally stable for <f> > ^j^lt^^i^ and l0Cally unstable otherwise while the lower endemic equilibrium becomes locally unstable for (f> > ^ ^ ^ + A ) ^ ^ 7 + ^ + / J ) A ^ • ^ N summary these equilibria exchange stability as the endemic equilibrium moves through the the disease-free equilibrium at <j> = ^ff^+A)^^7+a^+/x)A^ a n c ^ t n e r e exists only one epi-demiologically feasible endemic equilibrium for <f> < ^ ^ | ^ + A ) ^ ^ 7 + a ^ ) A ^ • ^ e e x i s t e n c e of a transcritical bifurcation at (/, </>) = (0, ^\il^^~+^A)]) using the theorem by So-tomayor (1973). Also these phenomena involving transcritical bifurcation is demonstrated in Figure 2.4. Theorem 2.3.2. (Sotomayor, 1973) Let x = f(x,p) (2.3.3) be the system of differential equations in W1 depending on the singular parameter p. When p = po, assume that there is an equilibrium XQ, i.e. f(xo,po) = 0. Suppose that the n x n matrix A = Df(xo,po) has a simple eigenvalue A = 0 with eigenvector v and that AT has an eigenvector w corresponding to the eigenvalue A = 0. Furthermore, suppose that A has k eigenvalues with negative real parts and (n — k — 1) eigenvalues with positive real parts (counting multiplicity) and that the following conditions are satisfied wTf^x0,po) ± 0, wT[D2f{x0,po)(v, v)] ± 0. (2.3.4) 39 Figure 2.4: Bifurcation of SIS model when there is no disease-fatality ,with various fraction infective from immigration, p. From the top left, various p values are used : a = 0,a = 0.05,a = 0.15 and a = 0.1 in clockwise direction respectively. Other parameters are fixed at p = 0, A = 3, 0 = 1, a = 0.2, p = 0.1, 7 = 12, 9 = 0.5 and A = 2. 40 Then there is a smooth curve of equilibrium points of 2.3.3 in Rn x Rn passing through (XQ,LIO) and tangent to the hyperplane Rn x {LIO}. Depending on the sign of the expressions in 2.3.4, there are no equilibrium points of 2.3.3 near Xo when LL < LLQ (or when LL > LIO) and there are two equilibrium points of 2.3.3 near XQ when LI > LIQ (or when LI < Lio)- The two equilibrium points of 2.3.3 near XQ are hyperbolic and have stable manifolds of dimension k and k + 1 respectively; i.e., the system 2.3.3 experiences a saddle-node bifurcation at the equilibrium point XQ as the parameter LL passes through the bifurcation value LL = LIQ. The set of C°°-vector fields satisfying the above condition is an open, dense subset in the Banach space of all C°°, one-parameter, vector fields with an equilibrium point XQ having a simple zero eigenvalue. Corollary 2.3.3. [12] If the conditions 2.3.4 a r e changed to wTffi{x0, Lio) = 0, wT[Dfli(x0,Li0)v} ^ 0 , and (2-3-5) wT{D2f(x0,Lio)(v,v))^0, then the system 2.3.3 experiences a transcritical bifurcation at the equilibrium point XQ as the parameter LI varies through the bifurcation value fj, = LLQ. 41 2.4 S p e c i a l cases 2.4.1 The case without disease fatalities in the absence of immigration, i.e. a = 0 and A = 0 Assuming that there is no immigration (A = 0) and that there are no disease fatalities. Let us investigate the new criterion which determine the presence of a backward bifurcation. We now can formulate the revised model: S' = A - 0SI - (LI + (f>)S + -yI + 6V I' = BSI + aBVI-{ii + i)I (2.4.1) V = <j>S- aBVI - (M + 0)V Since N = S + I + Vwe can add the three equations, which gives N' — S' + I' + V = A — LIN. 'By the theory of asymptotically autonomous systems, one can reduce the system to two dimensions by letting N = K and 5 = K — I — V: I' = 0[K - I - (1 - a)V]I -(n + 7 ) / (2.4.2) V' = <t>[K-I-V\- aBVI - (JJ, + 0)V One can solve for V using the first equation and substitute unto the second one to eliminate V, which yields EI2 + BI + C = 0 where E = a0 B = (n + 6 + a<j)) + <r(fi + 7) - a0K c=^ + ^ + 6 + ^-(, + e + a m Implicit differentiation with respect to 0 gives 42 Figure 2.5: Bifurcation of SIS model in the absence of immigration and when there is no disease fatality, i.e. a = 0, E = 0 = aI(K - /) + {p + 9 + a)(p + 9 + 4>) 02 The right hand side is greater than zero so we want 2EI + B > 0 to have the positive slope of bifurcation curve, and 2EI + B < 0 for a negative slope. Thus 2EI + B — 0 will be critical in determining the slope of the bifurcation curve. 43 2.4.2 The case when all susceptibles are vaccinated immediately and vac-cine does not lose its effectiveness i.e. 9 = 0, <fi —* oo Let us assume that there is no loss of effectiveness of vaccine , so we let 9 = 0. Also we assume that all the susceptibles are vaccinated immediately, i.e. <fi — • oo. If there is no fatality from the disease then now we can reduce a system to one-dimensional as below: /' .= PA + a/3I{K + i)I (2.4.3) = -a/312 + (a/3K - Li-j)I + pA At equilibrium, I' = 0 and let us assume that I*2 are roots of / ' = 0, with 1% > 0 and P2 < 0. „ _ - ( M + 7 - <T0K) ± y/{n + 7 ~ v0K)2 + 4a8pA h>2 ~ • 2*0 Let us define £ as £ = — (^  + 7 — CT0K). Now we can rewrite the endemic equilibrium as . l l ~ 2a~0 When p gets close to zero, one can study the limit of i | : l i m / ; = i ± f = ( ° i t { < ° ; 2"0 i f « > o Thus we can conclude that £ = 0 is a threshold that determines whether there exist an endemic equilibrium when p goes to zero. Since RQ = RQ = = GRQ < RQ. Note that RQ < RQ. Also if p goes to zero, I\ becomes zero, so there cannot be a backward bifurcation in this case. 44 The case with no infective immigrants In this subsection we make an additional assumption about immigrants. First, we assumed already that vaccine is applied immediately to susceptibles and that the vaccine does not lose its effectiveness. In addition to that, in this subsection, we also assume that there are no infective immigrants. Then we can let p = 0 in the equation 2.4.3, so the revised equation is /' = oQI{K -I)-{p + 7 ) I (2.4.4) = - a / 3 / 2 + {oBK -p,-i)I and its steady states are i ? = 0 /2- = ^ ; / " 7 (2.4.5) where I | is biologically relevant only when aBK — p — 7 > 0. It can be shown that equation 2.4.4 exhibits a transcritical bifurcation at /j* = 0 when the contact rate, 6, is used as a varying parameter. We use the theorem 2.4.1 to prove this. Theorem 2.4.1. Consider the function, x = f(x, p) , I , f l £ l . Let's assume that f has continuous derivatives up to third order for all (p,x) near the 45 bifurcation point (p*,x*). If there is a pair (p*, x*) for which f(x*,p*) = 0 fx(x*,p*) = 0 Mx*,p?) = 0 U(x*,P*) ? o fxx(x*,p*) + 0 then x — f(x,p) has a (two branch) transcritical bifurcation with quadratic tangency at (p*,x*), [12]. In our equation we let /(/, 0) = - a / 3 / 2 + {a(3K - p - ^ I Then fx(x, 0) = - 2 a 0 1 + (a0K - p - 7) fp(x,0) = -aI2 + aKI fx(){x,P) = -2aI + o-K fxx(x,0) = -2a 0. 46 0 1 2 3 4 5 0 1 2 3 4 5 Figure 2.6: Bifurcation when all susceptibles are vaccinated immediately and vaccine does not lose its effectiveness i.e. 9 = 0, <j> ~* 0 0 • Left: the case with no infective immigrant, i.e. p = 0 Right: the case with p = 0.1. Other parameter values are fixed at a = 0.2, p — 0.1, 7 = 12, A — 2, K = 20 while a control parameter is a contact rate, B. We can evaluate f(I,B) and these derivatives at (I*,0*) = (0, = o /*(**,/n = o /M**,/?*) = o fxp(x\B*) = o-K ^ 0 where other parameters are strictly positive. Hence according to the theorem 2.4.1 it follows that I = /(/,/3) has a (two branch) transcritical bifurcation with quadratic tangency at (I*,0*) = (0, ^ ) . This is demonstrated in Figure 2.6. 47 2.4.3 The case when the vaccine is utterly ineffective, i.e. a = 1, and there are no disease fatalities, i.e. a = 0 Let us assume that the vaccine is utterly ineffective, i.e. a = 1, and that there is no disease-related fatality, i.e. a = 0. One can modify the basic model, 2.1.1, based on this assumption in order to get the system: /' =pA + B(K - 1)1 + 7 ) I (2.4.6) V =</>(K-I)-BVI-(n + e + </>)I Note that this is somewhat similar to the previous system for the general SIS model. At equilibrium, we have -BI2 + {BK-Li--f)I+PA = 0 The possible real two roots of an equation are _ (BK-Li--y)±^/(BK-Li-j)2 + 4BPA h ' 2 ~ 2B From this expression, note that/I < 0 unless p=0 (so 4/3pa = 0) and II is positive if {BK-fjL-y) + v W - M - l)2 + 40pA > 0 1. If Q > 0 then I{ > 0 2. If Q < 0 then If < 0 when y/Q2 + AfipA > -Q i.e. AfipA > 0 where Q = /3K - /i - 7 In summary, if BK — LI — 7 > 0 then there always exists an endemic equilibrium in this case and no backward bifurcation is possible. 48 2.4.4 The case when the vaccine is absolutely effective a = 0 and no death is caused by the disease, i.e. a = 0 Let us assume that the vaccine is absolutely effective and no death is caused by the disease, i.e. a — 0. One can modify the basic model, 2.1.1, based on this assumption in order to get the system: I' =pA + 0[K -I - V]I - (/i + 7 ) / (2.4.7) v' = <t>[K -1] - {p + e + <j>)v We can get equilibrium conditions by letting the equation for / ' and V to be zero. pA + 0[K-I- V)I = (p + 7 ) / <t>[K-i\ = {p, + e + 4>)v From the second condition one can solve for V as following; V* 4>K-4>I P + O + 4> By substitution for V into the first condition gives PI*2 + QI* + R = 0 where 1. The case when p = 0 i j =0 r* — lo  P = Q = -0{p + 9 + 20) p + 9 + 4> 0K(p + 9) p + 9 + i R = pA - ( M + 7) -Q (3K(P + 0) p + 9 + (f> (p + l) P + 9 + , -0(p + 9 + 2(f>) Note that I\ > 0 if and only if 5 = BK(p + 0) - (p + j){p + 0 + 0) < 0 2. The case when p ^ 0 Equilibrium satisfy / ' = 0 and V = 0. M + 0 > + 0 + . > 0 with r _5±^82 + 4pA{p + 9 + 4>)f3{p + 0 + 20) 1 , 2 ~ 2/3(/x + 0 + 20) 5+\S\ P ™ 1 ~ 2/3(/x + 0 + 20) 0 if 5 < 0; if 6 > 0 ^ /3(p-+0+24>) Thus, 5 = 0 is a threshold, suggesting that there is a basic reproduction numbe RQ = 3K{p + 9) {p + 7 ) ( / J + 0 + 0) Note that lim /* = < p->0 0 if -Ro < 1; > 0 if Ro > 1 For p « 0 using Taylor's series (1 + x) 1 / 2 « 1 + §, 2A(JU + 9 + 0)/3(/J + 9 + 20) P l 2/3(^ + 0 + 20)/* « <5 + |5|[1 + <52 If RQ > 1 so that 5 > 0 this gives _ AQu + 0 + 0) 50 Thus at least for p close to zero, RQ represents a threshold. Figure 2.7: Bifurcation when vaccine is absolutely effective a — 0 and no death is caused by the disease, i.e. a = 0 Chapter 3 The SVIR model 3.1 T h e general S V I R M o d e l Let us now consider an SIR type disease when a vaccination program is in effect and there is a constant flow of incoming immigrants. We define S, V, I, R and N as the number of susceptibles, vaccinated, infectives and recovered individuals and total population respectively. We model new infections using the simple mass-action law, so that in general there are fiSI new infections in unit time when fi is the rate of contact that is sufficient to transmit the disease. We also assume a constant recovery rate, 7. The vaccine has the effect of reducing the susceptibility to infection by a factor of <j, so that a = 0 means that the vaccine is completely effective in preventing infection, while a = 1 means that the vaccine is utterly ineffective [10]. The rate at which the susceptible population is vaccinated is (f>, and the rate at which the vaccine wears off is 9. We assume that there can be disease related death and define a to be the rate of disease-related death, while fi is the rate of natural death that is not related to the disease. The population is replenished in two ways, birth and immigration. We assume that all newborns enter the susceptible class at the constant 51 52 rate of A and there is a constant incoming flow of immigrants, A, where some portion, p, of immigrants are infectives. We can now formulate this model, dividing the population into four classes - susceptibles (S), vaccinated (V), infectives(J) and recovered (R): S' = A + (1 - p)A - 0SI- (n + 0)5 + ev V' = <pS - afiVI - {n + 0)V (3.1.1) T =pA + BSI + aBVI - {LL + 7 + a)I R' = 7 / - LlR with S(t) + V(t) + I(t) + R(t) = N{t). Thus it follows that TV' = A + A - al - LiN. We can get an alternate but equivalent model by replacing SbyN — V — I — R. Now the model is V' = 4>(N - I - R ) - aBVI -{(p + e + Li)V I' =PA + 0I[N - I - ( l - a ) V - R ] - ( a + >y + Li)I (3.1.2) R' = 7 / - JJLR N' = A + K-al - LIN We can write the equilibrium conditions by letting the right-hand side of 3.1.2 to be Figure 3.1: The diagram for SVIR model S(t) Number of susceptibles at time t V(t) Number of vaccinated individuals at time t J(t) Number of infectives at time t R(t) Number of recovered people with immunity at time t A Immigrants P Portion of infectives among immigrants A Birth rate 0 Contact rate 7 Recovery rate <t> Vaccination rate a Factor by which the vaccine has the effect of reducing the infection e Rate at which the vaccine wears off p Natural death rate that is not related to the disease a Disease-related death rate Ro Basic reproductive number R* Vaccination reproductive number Table 3.1: Summary of notation 54 zero. The equilibrium conditions are • fiN = A + A - al pR = ^I [o0I + (cf> + 9 + p))V = <t>(N-I-R) (a + y + fj)I = pA + 0I[N - I - R - ( l - a)V) We can solve the first and the second conditions for the equilibrium values . .„ A + A - al* N = P il* M Solving for V from the third condition and substituting for N* and R* gives = ftN-I-R) (T0I(li + 9 + <f>) i/A+A-al _ j _ j£ a0I{ix + 9 + 4>) 4>[A + A - (a + 7 + u)I] n[o-pi + (n + e + <i>)] Eliminating TV, V, and R by substitution of these expressions into the last condition and simplifying, we obtain a cubic equation for the equilibrium values of / of the form EI3 + BI2 + CI + D = 0 (3.1.3) with E = 02a{a + 7 + p) B = a0p{a + j + p)- 02a(A + A) + 0(a + 7 + p){p + 9 + o(p) C = -pAap + p(a + 7 + p){p + 9 + 4>) - 0{p + 9 + a(p)(A + A) D = -pAp(p + 9 + a4>) 55 Let us first consider the basic reproductive number in the absence of vaccination. The basic reproductive number of the model without vaccination is * 0 = ^ ± f L (3.1.4) The reproductive number of the model with vaccination is y Y I fi(a + 7 + + 0 + <t>) 3.1.1 The vaccine reproductive number, R(<f>) Let us consider the case when there are no infective immigrants, i.e. p — 0. Now we consider the special case when p = 0 from the model 3.1.2. V' = 4>(N - I - R ) - apVI -((/> + 6 + u)V I1 = PSI + apVI - (a + 7 + n)I (3.1.6) R' = 7 / - fj,R N' = A + k - a l - LIN We can write the equilibrium conditions by letting the right-hand side of 3.1.6 be zero. The equilibrium conditions are o-pVI + ((j) + 6 + n)V = <t>(N - I - R ) pI[N - I - R - { I - a)V] = (a + 7 + n)I fj.R = jI LiN = A + A-al 56 We can solve the third and the fourth conditions for the equilibrium values R* ill (3.1.7) TV* = A + A-al* (3.1.8) From the second condition we can easily see that there is disease-free equilibrium, If = 0. At disease-free equilibrium we can evaluate the other equilibrium values of R and TV using 3.1.7 and 3.1.8 Rl = 0, J V * A + A A4 Solving for V from the first condition and substituting for TV* and R* gives V* 4>(N - I - R ) o-3I(n + 9 + <f>) i/A+K-al _ r ll G0I{LI + e + <t>) (j>[A + A - (a + 7 + fj,)I] H[<701 +(p + 0 + <l>)] Submitting the values of R, S and I at disease-free we find <t>{A + A) The Jacobian matrix of the system 3.1.6 is / -<f> - G0V •201 + 0[N - R - (1 - a)V] - (7 + a + u) J = V -a -o0I -{Li + e + (f>) -01(1 - a) 0 0 -4> 4>\ -01 01 -fi 0 0 -Ll/ 57 Now we can evaluate this Jacobian at the disease-free equilibrium -(Ll + 9 + cf>) -d> M 1 0 0 0 5 0 0 V —a 0 0 J1/r=0 The characteristic equation of the Jacobian matrix above is reduced to \0{A + A)(LI + 9 + O-(/>) fl(fl + 9 + (j>) (7 + a + fi) - A [-Gu+ (? + </> + A ) (/z + A ) 2 ] = 0 Thus, Ai A 2 A34 Ai(/i + 9 + <j>) -(li + 9 + ct>) - (7 + a + /i) Since A 2 and A 3 i4 are negative, the only condition required to ensure the stability of disease-free equilibrium is 0{A + A)(fi + 9 + ac/)) fi(fi + 9 + cj)) Thus now we define a vaccine reproductive number (j + a + fi) < 0 R(</>) = 0{A + A)(fi + 6 + acj>) fi('y + a + fi)(fi + e + 4>) and the disease-free equilibrium is locally asymptotically stable if R(<f>) < 1. In the absence of vaccination, i.e. </> = 0, we define a dimensionless basic reproductive number, RQ as 0(A + A) fi{j + a + fi) 58 Note that R(0) = Ro and that R(cf)) < Ro and also that aRo < R((f>) for all cf> > 0 since 0 < a < 1 with lim R((f>) = aRo-4>—>oo 3.2 E n d e m i c e q u i l i b r i a a n d b i f u r c a t i o n s w h e n t h e r e a re n o i n f e c t i v e i m m i g r a n t s , p = 0. Let us first consider the two extreme cases in order to investigate endemic equilibria for the model 3.1.2. The first case is when the vaccine is useless and there is no infective immigrant, i.e. a = 1 and p — 0. Our system becomes S' = A + A - BSI - {p + <f>)S + 9V V = <f>S- BVI -{p + 9)V I1 = BSI + 8VI - (M + 7 + <*)! RJ = jl — pR with S(t) + V(t) + I(t) + R(t) = N(t). Note that R((j>) reduces to R0 when a = 1. By substitution and using'the equilibrium conditions, one can find there is an endemic equilibrium, r = B(A + A)- pja + j + p) Pia + j + p) which exists only when RQ > 1. If instead we suppose that the vaccine is completely effective and that there are no 59 infective immigrants, i.e. a = 0 and p — 0, then the model becomes S' = A + A - 0SI - {p + <f>)S + 9V V = 4>S-{p + 9)V I' = 0SI-{n + i + a)I R' = 7 / - ^i? with S(t) + V(t)+I(t)+R(t) = N(t). After complicated calculation which we omit we find that there is one endemic equilibrium, r = 0(A + A){0 + p)-p((f> + 9 + p){a + p + 7) 0{a + p, + i){9 + p) which exists only for R{4>) > 1 when a = 0. Now we wish to consider the more general case when the vaccine is partially ineffective and when there are no infectives immigrant, i.e. 0 < a < 1 andp = 0. S' = A + A - 0SI - (p + <t>)S + 9V V' = <f>S- (70VI -{p + 9)V (3.2.1) I' = 0SI + <70VI -(p + j + a)I R' = 7 / - u.R with S(t) + V(t)+I(t)+R(t) = N(t). From equilibrium condition 3.1.3, we know that there is a disease-free equilibrium regardless of different parameter values. After factoring out this disease-free equilibrium, / , we get an equilibrium condition as a quadratic polynomial of/. EI2 + BI + C = 0 (3.2.2) 60 with E = 02a{a + 7 + M) B = o0fi(a + <y + fi)- 02a(A + A) + 0(a + 7 + fi)(fi + 9 + acp) C = fi{a + 7 + u)(n + 9 + (f>) - 0(fi + 9 + a<p){A + A) Now an endemic equilibrium corresponds to a positive real solution of equation 3.2.2. Note that E > 0 and that C < 0 precisely when R(cf>) > 1. Note also that B2 - 4EC > 0 when C < 0. One can easily deduce that there is precisely one endemic equilibrium when R(4>) > 1 since there are two real roots and the product of those two roots is negative. On the other hand we can see that C > 0 if R(4>) < 1. Note that there are exactly two changes in the sign of coefficients of equation 3.2.2 if coefficient B < 0 and none when B > 0. By Descartes' rule of signs one can conclude that the maximum number of endemic equilibria is two when R(<f>) < 1 and B < 0, and that there is no endemic equilibrium when R((f>) < 1 and B > 0. However it is shown that it is always the case that the system does not have any endemic equilibrium when R((/>) < 1. Proposition 3.2.1. . For model 3.2.1, with R(<j>) as defined in 3.1.5, there is no endemic equilibrium when R(<p) < 1. Proof. We first assume that R(<p) < 1 and B < 0: R(<j>)<l & 3(n + 9 + a(f>)(A + A)<iJ,(a + <y + iJ.){iJ, + 9 + (l>) B<0 & 0o~fx(a + 7 + /i) + 0(a + 7 + n){u- + 9 + a<p) < 32a(A + A) Combining these two conditions one can get the following relation: o\i(a + 7 + u)(u- + 9 + cf>) a[i{a + 7 + n) + (a + 7 + u)(u- + 9 + a<p) < fi + 9 + aq 61 Factoring out a + 7 + LI it reduces to ULL + 11 + 6 + a4> < Ll + 9 + a<p Multiplying LI + 6 + <J(f> on both sides of equation gives [(a + 1)LI + 9 + O-4>]{LI + 6 + a<t>) < GLL^LI + 9 + $) After expansion and some calculations which we omit we get an expression, (9 + (7<J>)(LI + 9 + acp) + LL{LL + 9) + a 2 ^ < 0 which is a contradiction for all nonzero parameters. Therefore we can rule out the case where R{4>) < 1 and B < 0. Note that R(<f>) < 1 corresponds to C > 0 where B and C are defined as in 3.2.2. One can deduce that when R(<f>) < 1 73 is always positive for nonzero parameters, thus we cannot have any endemic equilibrium in this case by Descartes' rule of signs. • Theorem 3.2.2. Descartes' rule of signs The number of positive real roots of an equation with real coefficients, f(x) = a0xn + aixn~x -I h an = 0, is never greater than the number of variations in the sequences of its coefficients ao, 01, • • • , On and, if less, then always by an even number We can now summarize the results as follows: Lemma 3.2.3. For model 3.2.1, with R(4>) as defined in 3.1.5, 1. When R((f>) > 1, there is precisely one endemic equilibrium. 62 2. When R(4>) < 1, there is no endemic equilibrium. 3.3 The general SVIR case: Non-existence of backward bi-furcation Now we return to more general case of the SVIR model, 3.1.1 and we wish to study on the types of its bifurcation. We recall the endemic equilibrium condition as a cubic equation of I: EI3 + BI2 + CI + D = 0 with E = 02o(a + 7 + LI) B = o0Li(a + 7 + LI) - 02o-{k + A) + 0{a + 7 + u)(u- + 9 + acp) C = -pAo-Li + Li(a + 7 + + 9 + (f>) - 0(LI + 9 + <x<p)(A + A) D = -PALI(LL + 9 + o-<j>) Note that it is clear that E > 0 and D < 0. As the vaccination rate <j> —>• 00, the reproductive number approaches CT.RO- A backward bifurcation is signalled by an interval for cf> on which there are three positive roots of 3.1.3. In order to have a backward bifurcation at R(<j>) = 1, it is necessary to have R0 > 1 aRo < 1 R{<f>) < 1 63 Under the assumption aRo < 1 we have B = a0n(a + u- + 7) - 82a{A + A) + 0(a + p + 7 ) ^ + 9 + acj>) > a0n(a + \i + 7) - 0u-{a + \i + 7) + 0{a + \i + j)(fj, + 9 + a<p) = 0(a + (j, + 7)(cr/j + 9 + a<f>) > 0 Thus regardless of the sign of C there is only one sign change in the coefficients of 3.1.3. By Descartes' rule of signs the equation 3.1.3 can have at most one positive root. Thus there cannot be a backward bifurcation. 3.4 O u t g o i n g immigran t s and backward bifurcat ions In previous section we showed that there cannot be a backward bifurcation for all nonneg-ative parameter values including immigration rate. However, from numerical simulation using XPPAUT, it was found that there can be backward bifurcation if we allow the rate of immigration to be negative, i.e. A < 0. In other words, if there is a constant outgoing immigration flow from the host population, backward bifurcation can occur with some pa-rameter sets. Some exemplary figures are shown below with indication of Hopf bifurcation and its stability. For instance, in Figure 3.3, Hopf bifurcation occurs at <f> = 0.4727 with period 4.557 and bifurcation loses stability at that point to become unstable as <fi increases through 0 = 0.4727. 64 Figure 3.2: Outgoing immigrants and Hopf bifurcation 1: Other parameters are fixed at A = -5.55, p = 0.50, a = 0.43, 0 = 2.92, a = 0.34, u- = 0.68, 7 = 0.53, 9 = 0.36 65 Figure 3.3: Outgoing immigrants and Hopf bifurcation 2: Other parameters are fixed at A = -7.65, p = 0.10, a = 0.73, 0 = 3.77, a = 0.18, p = 0.37, 7 = 0.92, 6 = 0.52 66 Figure 3.4: Outgoing immigrants and Hopf bifurcation 3: Other parameters are fixed at A = -2.73, p = 0.47, a = 0.90, 0 = 6.26, a = 0.06, LI = 0.09, 7 = 2.71, 0 = 0.41 67 3.5 Spec ia l case: T h e S V I R m o d e l w h e n there are no disease fatalit ies but w i t h infective immigran t s , i.e. p / 0 and a = 0 In this section we assume that there is a constant flow of infective immigrants into host population and no disease fatality. Based on this assumption, we have the system of ODE's below: S' = A + (1 - p)A - 0SI - (LI + <t>)S + 9V l'=pA + BSI + aBVI -(fj. + 7 )7 V = (t>S-aBVI-(Li + 0)V R' = 7 / - / i i i The total population is the sum of population of each class, i.e. N = S + I + V + R. Thus N' = S' + I' + V + R = A + A — fiN. By theory of autonomous system, it follows that limt_oo N = K — Replacing S with K — I — V — R,we have reduced system of ODE's here: • l'=pA + fi{K - I -V - R)I + aBVI -(fi + 7 ) / V' = 4>(K - I - V — R) - aBVI -{9 + LL)V (3.5.1) R' = 7 / - fiR 68 Equilibrium Condition Endemic equilibrium conditions are ^ - 0 1 = -O-0V + {p + 7) - 0{K - V - R ) 4>(K - I — V - R) = O-0VI + {9 + p)V 7 / = pR One can reduce these endemic equilibrium conditions into one cubic equation of / by sub-stituting for K,V and R: f(I) = EI3 + BI2 + CI + D = 0 (3.5.2) where E = B2a{p + 1) B = 0[-0a(A + A) + (p + 7)(0 + P + o<t> + ap)} C = {-pApo-0 - B(A + A)(9 + p + acf>) + p{p + + P + <f>)] D = —pAp2(9 + p + <p) Since 3.5.2 /(0) < 0 and l im^oo = 00 there exist either one or three positive roots, I*. Now let us consider /'(/) = 3EI2 + 2BI + C. By Rolle's theorem, if /(/) has three distinctive positive rots then /'(/) must have two positive roots. Thus B < 0 and C > 0 is a necessary condition to three positive endemic equilibria. 69 3.5.1 Stability analysis Linearization of 3.5.1 will lead to / 0{K - V - R ) - 2 0 I + apV - ( M + 7) -01 + a0I -01 \ -4> - (70V —(f) — 9 — Ll — (701 -(f) \ 7 0 -LI J (3.5.3) 3.5.3 can be reduced using equilibrium condition, ^+01 = 0{K — I — R) +o0V — (/z + 7) , to the following: (-P-T~0I 0(a - 1)1 -0l\ IA B c\ -(j) - (70V -{(f> + 9 + Ll)-(70I -<p = D E F \ 7 0 -LL) \G 0 The characteristic equation is A 3 + aiA 2 + a 2A + a 3 = 0 where ai = -{A + E + I) a2 = [EI + AE + A I - C G - BD) a 3 = {BDI + CGE - AEI - BFG) For global stability, it is required that ai > 0 a 3 > 0 a\a2 — a 3 > 0 Chapter 4 Conclusion and Discussion The purpose of this paper is to take a close look at the endemic behavior of the diseases of SIS and SIR type. To a simple SIS model with vaccination we added the immigration of infectives and the rate of the disease-related death rate. As to the contact between infectives and susceptibles we assume a bilinear incidence. The result of out mathematical analysis indicate that a vaccination campaign <f> has an effect of reducing a reproductive number, which means that the average number of secondary infection caused by an average infective becomes smaller when vaccination is in effect. Furthermore, in SVIS model, a vaccination campaign meant to reduce a disease's reproductive number below one may fail to control disease when there is a backward bifurcation. It one aim to eradicate the disease that is already established, he needs to make sure the vaccination rate is far greater than the one corresponding the vaccination reproductive number, R((p) — 1. Bringing down the vaccina-tion reproductive number just below one may not be good enough to eradicate the disease in such a case. Also if there is no immigration of infective, a typical transcritical bifurcation may be observed. The disease-free equilibrium and endemic one coincide at R(<f>) = 1 and they exchange the stability at that point. Also we found that there is possible forward 70 71 bifurcation with exchange of stability using some parameter sets. The main results from the SVIS model is that a backward bifurcation is possible if the vaccination is partially but not completely effective and for the SVIR model backward bi-furcation is not possible. Note that both SVIS and SVIR model are set up based on bilinear incidence and that standard and more general incidence is not covered, but if there are no disease deaths so that the total population size is (asymptotically) constant, all incidence forms are equivalent. From SVIR model we proved that there cannot be a backward bifurcation and only for-ward one is possible with incoming infective immigrants. However in very rare cases there is possible backward backward bifurcation only with outgoing immigration of infectives. The model can be more significant when transformed to adjust to specific transmittable diseases; using specific parameter values and adding some extra terms into a model would do this. If we can get a real clinical data then we will be even able to compare one to the results. Furthermore exposed periods of infection can be also considered to improve a model. Regarding infective periods, I can improve it by making it non-fixed but arbitrary; some studied on exponentially distributed one, for example. Also the contact rate with saturation may be used in some cases. 72 Speaking of tools for this study, I mainly dealt with a system of nonlinear ordinary dif-ferential equations and understand it by looking at equilibrium, phase plane and bifurcation diagram. Stability analysis is critical in this study since we would be able to know whether endemic equilibrium would be stable so the disease would persist or not. By setting up a good epidemic model and understanding well we can have many advantages of preventing invasion of infection to population. One of the most important result may be knowing the possible endemic equilibrium. To prevent it, it will certainly be helpful for us to plan vac-cination policy better. Knowing bifurcation point of backward bifurcation we would know the minimum portion of susceptible class to be vaccinated in order to have a control over a disease. Also having sufficient condition for having backward bifurcation would help us to avoid sudden jump in number of infectives due to unstable branch of bifurcation [10]. Appendix A Bifurcation A . l H o p f B i f u r c a t i o n Theorem A.1.1. (Hopf [1942]).Suppose that the system x = f^x), x G I K " , p £R has an equilibrium (xo,po) at which the following properties are satisfied: (HI) Df(xo) has a simple pair of pure imaginary eigenvalues and no Then (HI) implies that there is a smooth curve of equilibria (x(p),p) with x(po) = XQ. The eigenvalues X(p), X(p) of Dxftl0(x(p)) which are imaginary at p = po vary smoothly with p. If, moreover, A ( i ? e A ( M ) ) | M = M 0 = d^0, then there is a unique three-dimensional center manifold passing through (XQ, po) m 1" x R and a smooth system of coordinates (preserving the planes p = const.) for which the Taylor expansion of degree 3 on the center manifold is given by x = (dp + a(x2 + y2))x - (to + cp + b(x2 + y2))y, (A. l . l ) y = [u> + cp + b{x2 + y2))x + (dp + a(x2 + y2))y. If a ^ 0 there is a surface of periodic solutions in the center manifold which has quadratic tangency with the eigenspace of \{po), MPO) agreeing to second order with the paraboloid p = —(a/d)(x2 + y2). If a < 0, then these periodic solutions are stable limit cycles, while if a > 0, the periodic solutions are repelling. ' Theorem A.1.2. Suppose that system with x G W1 and p G W1 has a critical point XQ for p = po and that Df(xo,po) has a simple pair of pure imaginary eigenvalues with zero 73 74 real part. Then there is a smooth curve of equilibrium points x(fi) with X(LIQ) = XQ and the eigenvalues, \(fi) and X(fi) of Df{x(LI), fi), which are pure imaginary at LL = LLO, vary smoothly with LI. Furthermore, if then there is a unique two-dimensional center manifold through the point (xo,(io) and a smooth transformation of coordinates such that the system on the center manifold is trans-formed into the normal form x = -y + ax(x2 + y2) - by(x2 + y2) + 0(\x\4) y = x + bx(x2 + y2) + ay{x2 + y2) + 0(|x| 4) in a neighborhood of the origin which, for a ^ O , has a weak focus of multiplicity one at the origin and x = fix — y + ax(x2 + y2) — by(x2 + y2) y = x + fiy + bx(x2 + y2) + ay(x2 + y2) is a universal unfolding of this normal form in a neighborhood of the origin on the center manifold. Here we summarize theorem for some types of bifurcations, saddle-node, transcritical, and pitchfork bifurcation of x = f(x, fi) , i , / i £ l T h e o r e m A . 1 . 3 . Consider the function, x = f(x, fi) , x,fiGR. Let's assume that f has continuous derivatives up to third order for all (fi,x) near the 75 bifurcation point (p*,x*). If there is a pair (p*,x*) for which f(x*,p*) = 0 fx(x*,p*) = 0 U(x*,p*) ± 0 fxx(x*,P*) ? 0 then x = f(x,p) has a saddle-node bifurcation with quadratic tangency at (p*,x*). T h e o r e m A.1 . 4 . Consider the function, x = f(x, p) , x,p €M.. Let's assume that f has continuous derivatives up to third order for all (p,x) near the bifurcation point (p*,x*). If there is a pair (p*,x*) for which f(x*,p*) = 0 fx(x*,p*) = 0 = o U(x*,P*) ^ o fxx(x*,p*) ^ 0 then x = f(x,p) has a (two branch) transcritical bifurcation with quadratic tangency at (p*,x*). T h e o r e m A . 1 . 5 . Consider the function, x = f(x, p) , x,p GR. Let's assume that f has continuous derivatives up to fourth order for all (p,x) near the bifurcation point (LI*,X*). If there is a pair (LL*, X*) for which f(x*,[L*) = 0 fx(x*,(l*) = 0 U(X*,LL*) = 0 fxx(x*,U*) = 0 fxfi(x*,Ll*) # 0 fxxx(x*,Ll*) ± 0 then x = f(x,Li) has a pitchfork bifurcation with quadratic tangency at (LI*, Appendix B Transcritical bifurcation of SVIS model with p = 0 Our system of SVIS type with no infective immigrant is /' = 0I[N - J — (1 — a)V] - [fi + 7 + a) I V' = <j>(N - I ) - o-BVI -{fi + 9 + (p)V N' = A + A- fiN - a l with N{t) = S(t) + I(t) + V(t) and N'(t) = S'{t) + I'{t) + V'{t) = A + A - fiN{t) - al(t). In the notation of Theorem by Sotomayor, we have A = Df(I — 0,(f> = 0o) / 0 0 0 \ V -a (B.0.1) 77 78 where = (LI + 9)[LI{J + a + LI) — 0(A + A)] 9 0 o0(A + A) - ii(7 + a + LI) ( / 0(7 = O,0 = ^o) = -o-/3(A+A)+^(7+Q+M) V 0 V = Li + e + fo w = (l,0,0)T {LL + a)(ii + 9 + M + a0(A + A) fl(jM + 9 + </>0) -)' (B.0.2) Appendix C Derivation of Stability Formula of two-dimensional Hopf bifurcation [5] If the reduced system has a purely imaginary pair of eigenvalues A = ±iu>, then we can represent it as a single complex expression, z = \z + h(z,z), (C.0.1) where z = x + iy and A = iu At [x = 0, the normal form, x = (d/z + a(x2 + y2))x - (u + CLI + b(x2 + y2))y, • y = (to + cu + b(x2 + y2))x + (dfj, + a(x2 + y2))y 79 80 becomes w = Xw + ciw2w + c2w3w2 -\ t-ckwk+1wk+ 0(\w\2k+3) (C.0.2) = Xw + h(w,w), (C.0.3) where the coefficients are of them form Cj = a,j + ibj. Note that in polar coordinates we have r = o i r 3 + a2r5 H , 6 = to + bxr2 + b2rA -| , and the first nonzero coefficients CLJ , bj determine the stability of the periodic orbit and the amplitude dependent modification to its period. We use the near identity transformation, z = w + ip(w,w), tp = 0{\w\2), (C.0.4) , in order to transform C.0.1 into C.0.2. Substituting C.0.4 into C.0.1 and using C.0.2, we obtain X(wtpw — tp) + Xwtpu, — h(w + tp, w + tp) — h(w, w)(l + tpw) — h(w, w)ipw- (CO.5) We proceed by expressing tp as a Taylor series tp{w,w)= ^fc^f+0(H4)- (C.0.6) 2<j+k<3 where ^•fc = dtfP+k/dwjdwk. 81 Now let's use the normal form, h{w,w) = ciw2w + 0(\w\5), and substitute CO.6 into CO.5 in order to get 9 -9 9 - 9 U! - _ - U) W _ VJ ^Pww-^- + Xipwwww + (2A - \)4>ww— = hww—+ hwidwvj + hyjW — +0(\w\3). Equating coefficients yields the leading terms in the transformation , h w w ihww . . ipww = —r- = , (C.0.7) A LO tpww = = , (C.0.8) A to fe=(2TrT) = ( a o - 9 ) We do expansion to one higher order and equate the coefficients of the normal form term uP"w. For this term, the coefficient on the left-hand side of CO.5 vanishes and the right-hand side becomes or, using C.0.7 ci = —(hwwhWw — 2|/Yu,ui|2 — -\hww^) H . (CO.10) Thus it follows that 2ai = 2Re(c{) = h*wW - - ( h ^ h 1 ^ + h ^ h ' j , (C.0.11) where R and I denote real and imaginary parts, respectively. Expanding / and g in Taylor series and taking real and imaginary parts of the complex 82 valued function h it follows that the relevant terms in CO.11 may be expressed as hwww = g ( / x x x + fxyy + 9xxy + 9yyy\ h-ww — ~^{fxx ~ fyy + ^9xy), hWw — ^{9xx 9yy 2fxy), hww = T(/X X + fyy). (C.0.12) h-ww — ^ (dxx + 9yy)-We can derive the stability formula CO.13 by substitution of C.0.12 into CO.11: 16o.i — ( / x x x " f " fxyy ~\~ 9xxy ~T" 1 "f~ [ / x y ( / x x ~t~ fyy) 9xy(^9xx ~f" 9yy) fxx9xx ~f" fyy9yy\-LO or, O l — [ / x x x ~r~ / x y y ~r" 9xxy ~\~ 9yyy] ~r~ [fxyifxx "t" / y y ) 9xy{,9xx ~\~ 9yy) fxx9xx ~t~ fyy9yy\-1UL0 (C.0.13) Also note that the normal form for the parameterized Hopf bifurcation is neatly expressed in complex variables as w = Xw + ciw2w + 0{\wf), (C.0.14) where A = LL + iu. Bibliography [1] F. Brauer and C. Castillo-Chavez 2001 Mathematical Models in Population Biology and Epidemiology. Springer Verlag, New York. [2] F. Brauer and J.A. Nohel 1989 The qualitative theory of ordinary differential equations: an introduction. Dover Publications, New York. [3] F. Brauer and P. van den Driessche 2001 Models for transmission of disease with immigration of infectives, Math. Biosci, 171(2):143-1542. [4] C. Castillo-Chavez and H. Thieme Asymptotically Autonomous Epidemic Models, Mathematical Population Dynamics: Analysis of Heterogeneity Vol. One: Theory of Epidemics (0. Arino, D. Axelrod, M. Kimmel, M. Langlais; eds.), 33-50. Wuerz 1995 [5] J. Guckenheimer and P. J. Holmes Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields, Applied Mathematical Sciences, Vol 42, Springer-Verlag, New York, 1983 [6] H.W. Hethcote 1976 Qualitative analyses of communicable disease models, Math. Biosci., 28(3/4):335-356. [7] E. Hopf 1942 Abzweigung einer periodischen Losung von einer stationaren Losung eines differential-systems. Ber. Math.-Phys. Kl. Sachs Acad. Wiss. Leipzig, 95(1):3-22. (A translation of this paper, and comments on it appear as Section 5 of Marsden-McCracken [1976].) 83 84 [8] D.W. Jordan and P. Smith 1977 Nonlinear ordinary differential equations. Oxford University Press, Oxford. [9] W.O. Kermack and A . G . McKendrick (1927) A contribution to the mathematical the-ory of epidemics, Proc. Royal Soc. London, 115:700-721. [10] C M . Kribs-Zaleta and Z.X. Vekasco-Hernandez 2000 A simple vaccination model with multiple endemic states Mathematical Biosciences 164: 183-201 [11] J.D. Murray 1989 Mathematical Biology, Biomathematics 19, Springer-Verlag, Berlin-Heidelberg-New York. [12] L. Perko Differential equations and dynamical systems, Texts in Applied Mathematics, Vol.7, Springer-Verlag, New York, 1983. [13] H.E. Soper 1929 Interpretation of periodicity in disease prevalence, J. roy. Stat. Soc, Series B, 92:34-73 [14] J. Sotomayor 1973 Generic bifurcations of dynamical systems. In Dynamical Systems, M . M . Peixoto (ed.), pp. 549-560. Academic Press: New York. [15] S. Strogatz 1994 Nonlinear Dynamics and Chaos, Addison-Wesley, Reading, Mass. [16] S. Wiggins, Introduction to Applied Nonlinear Dynamical Systems and Chaos, Texts in Applied Mathematics, Vol.2, Springer-Verlag, New York, 1990.