Theoretical epidemiology and parasitology Luděk Berec Department of Theoretical Ecology, Institute of Entomology, Biology Centre ASCR Institute of Mathematics and Biomathematics, Faculty of Science, University of South Bohemia Branišovská 31, 37005 České Budějovice, Czech Republic [email protected] www.entu.cas.cz/berec

Source: www.dshs.state.tx.us/idcu/disease/influenza/surveillance/2008/

Source: Keeling and Rohani (2008)

Source: www.aps.anl.gov/Science/Highlights/2005/20050608.htm

Source: jwit.webinstituteforteachers.org/∼naphhoff/htc213/Outline.htm

Aim of the course

To provide an introduction to theoretical epidemiology and parasitology, including development of mathematical models for the spread and control of infectious diseases, tools for their analysis, as well as their application potential

Literature

Brauer et al. (2008) Mathematical epidemiology

Keeling and Rohani (2008) Modelng infectious diseases

Anderson and May (1992) Infectious diseases of humans

plus many scientific papers

Introduction

Agenda for this part Mathematical modeling What is parasite? Why modeling of infectious diseases? Basic concepts in epidemiology Basic reproduction number R0 Outline of the course

Mathematical modeling

In science, model is a representation of (a part of) reality that we aim to study Physical models: material, pictorial, or analogical representations of (a part of) an actual system, such as scale material models used in wind tunnel experiments or flight simulators Mathematical models: actual systems are represented by mathematical objects and their relationships, often in the form of various types of equations or governing rules run as computer simulations

Growth of a population How a population changes in size with time? Mathematical representation: Time t, population size N(t) b . . . number of births per year per individual d . . . number of deaths per year per individual Dynamics: Year t: N(t) Year t + 1: N(t + 1) = N(t) + bN(t) − dN(t) = (1 + b − d)N(t) Year t + 2: N(t + 2) = N(t + 1) + bN(t + 1) − dN(t + 1) = (1 + b − d)N(t + 1) = (1 + b − d)2 N(t) N(t + k) = (1 + b − d)k N(t) b > d . . . population grows b < d . . . population declines

Mathematical modeling

Mathematical models are useful experimental tools for Building and testing theories Generating hypotheses Answering specific questions Determining system sensitivity to parameter changes Estimating key parameters from data Supplanting experiments that we, for some reason, cannot conduct practically Assessing and comparing various management actions before they are actually employed

Mathematical modeling Basic steps in developing mathematical models 1

Formulate the question

2

Determine the basic ingredients

3

Qualitatively describe the examined system

4

Quantitatively describe the examined system

5

Solve the model

6

Perform checks and balances: if not passed, refine the model and repeat steps 2-5

7

Relate the results back to the question

Contribution of the original scientific discipline, such as biology, can mostly be detached from contribution of mathematics

Mathematical modeling Any model is necessarily a simplified representation of reality We consider only those aspects of an actual system that are essential for its understanding and/or prediction of its future behavior No model is the best one, there is always a place for improvement, no matter how large From no model we can require more than conditional predictions of the type “what effects would a given situation imply if it occurs” Models just reveal consequences of the assumptions upon which they are based

Mathematical modeling

Mathematical modeling is always a trade-off between Simple (generic, strategic) models: highlight qualitative behavior Detailed (specific, tactic) models: specific quantitative predictions One should always start simply, and gradually add complexity so long as complexity also adds insight A model should be as simple as possible but no simpler

What is parasite?

Epidemiology deals with parasites Parasites spend most of their life in a host-parasite association, or alternatively only a short periods of time, adopting a free-living mode for the major part of their life cycle To be considered parasitic, a species should satisfy three conditions: Utilize its host as a habitat Nutritionally depend on its host Cause harm to its host

What is parasite?

Variability of parasites in the degree of harm to the hosts: Lethal infections – host deaths occur due to parasite infection, these will also kill the parasites contained within Symbionts – cause negligible, if any, harm to the host even if present in very large numbers Anything in between these two extremes

What is parasite?

Effects of parasites on the hosts vary: Increased mortality rate of the affected hosts (lethal diseases) Reduced birth rate, which can even reach the state of parasitic castration Infected individuals with lowered reproduction might live longer

Induced morphological changes in hosts, induced changes in host behavior, and biases in the host’s sex ratio at birth – all serve to increase transmission efficiency of parasites

Why modeling of infectious diseases? General interest in infectious diseases Identification of new diseases (Lyme disease, HIV/AIDS, SARS, swine flu) Frequent reappearance of existing diseases (plague, cholera) Emergence of antibiotic-resistant strains (tuberculosis) Spread into new regions due to climate changes (malaria) Opportunities for new and existing infectious diseases due to human or animal invasions of new ecosystems, increased international travel, and changes in social and economic patterns Recent popular books, movies and TV series have given us exciting accounts of the emergence and detection of new diseases

Why modeling of infectious diseases?

Experiments in epidemiology often difficult or impossible to design Serious ethical questions involved in withholding treatment from a control group Collected data often incomplete or inaccurate Mathematical models complement empirical work and are extremely important in: Understanding fundamental mechanisms that drive disease spread Suggesting strategies for disease control

Why modeling of infectious diseases?

More specifically, mathematical models of infectious diseases Help clarify assumptions, variables, and parameters Pathways involved in parasite spreading, degree of heterogeneity needed

Provide conceptual results such as thresholds for disease invasion and plausibility of parasite eradication Contribute to design and analysis of epidemiological surveys, especially by suggesting crucial data that should be collected Contact tracing for sexually transmitted diseases

Why modeling of infectious diseases?

More specifically, mathematical models of infectious diseases Can be used as experimental tools for testing control measures and determining sensitivities to parameter changes Design of vaccination strategies

Can be used to compare and optimize costs and efficiency of various detection, prevention and control programs When parasites are used as biocontrol agents, models provide insights into the circumstances under which parasites are capable of regulating a host population, and of doing so in a required and stable manner Myxomatosis to control rabbits in Australia

Why modeling of infectious diseases? Long history of disease modeling in mathematical biology Sir Ronald Ross ∼ 1900-1910’s, malaria

W. O. Kermack and A. G. McKendrick ∼ 1920-1930’s Recently, modeling has become a part of epidemiology policy decision-making in several countries (UK, USA, Canada) Predictions of models of HIV/AIDS, BSE, foot and mouth disease, measles and SARS have had an impact on public health policy in these countries A tremendous variety of models have already been developed, analyzed, and applied to a tremendous variety of infections, such as malaria, rabies, syphilis, Lyme disease, and virtually any childhood disease

Why modeling of infectious diseases?

Source: Keeling and Rohani (2008)

Why modeling of infectious diseases?

We will never be able to predict the precise course of a disease, or which individuals will be infected. The best that we can hope for is models that provide confidence intervals on the disease behavior and determine the risk of infection for various groups of hosts.

Basic concepts in epidemiology

Two broad groups of parasites Microparasites viruses, bacteria, fungi, some (parasitic) protozoa Macroparasites some (parasitic) protozoa, (parasitic) helminths and arthropods Differences in their life cycles and impacts on their hosts call for development of structurally different epidemiological models

Basic concepts in epidemiology

Microparasites Small size – generally cannot be seen with the naked eye Short generation times Extremely high rates of reproduction within hosts Tendency to induce immunity to reinfection Duration of infection typically short and transient Can complete full life cycle inside a single host

Basic concepts in epidemiology

Macroparasites Large enough to be seen with the naked eye Much longer generation times Multiplication in hosts either absent or occurs at a low rate Immune responses tend to be of a relatively short duration Macroparasitic infections tend to be persistent, with hosts being continually reinfected Typically need two or more host species to complete full life cycle

Basic concepts in epidemiology

To complete full life cycle, parasites may pass from one host to the next either directly or indirectly via one or more intermediate host species Intermediate hosts are hosts in which the parasites reproduce asexually or in which the parasitic larvae grow; mostly invertebrates, but also vertebrates, including humans Definitive hosts are hosts in which the parasites mature and sexually reproduce (if they have the ability to do so) and in which the parasite’s life cycle begins as well as ends; a definitive host can be the only host of a parasite

Basic concepts in epidemiology

Direct transmission Physical contact between hosts (sexually transmitted diseases) Transmission stages of a parasite leaving a host and then being picked up by another host through: Inhalation (influenza) Ingestion (pinworm) Penetration of the skin (hookworm)

All these paths are examples of horizontal transmission Vertical transmission: infection passes from a parent to its unborn offspring through egg, embryo or host chromosomes (e.g. HIV/AIDS and many viral infections of arthropods)

Basic concepts in epidemiology

Indirect transmission Occurs when life cycle of a parasite involves one or more intermediate hosts Biting by vectors (mosquitoes, flies, ticks) that serve as intermediate hosts Penetration by free-living transmission stages that are produced by other (e.g. molluscan) intermediate hosts Parasites can be ingested when an infected intermediate host is eaten by a predatory or scavenging definitive host

A parasite life cycle

Source: cs.wikipedia.org/wiki/Soubor:Schistosomiasis Life Cycle.jpeg

Basic concepts in epidemiology

Infectious diseases can display two different temporal patterns: epidemic and endemic Epidemic diseases Act on a short time scale Sudden outbreak of a disease that infects a substantial portion of the population in a region before it disappears 2002 epidemics of SARS, Ebola and bird flu Usually leave many individuals untouched Epidemics often recur, possibly diminishing in severity as populations develop some immunity (seasonal flu)

Basic concepts in epidemiology Endemic diseases Disease is always present Relatively small fluctuations in monthly cases counts Only slow increase or decrease over the course of years Malaria, typhus, cholera, and many other diseases are endemic in many parts of the world Intermediate scenario: diseases constantly present but outbreaking frequently Stochastic effects can play a role in classifying a disease as endemic or epidemic – endemic diseases the prevalence of which is so small as to give a high probability of stochastic fade-out can actually be viewed as recurrent epidemics

Basic reproduction number Quantity of central importance in epidemiology Traditionally denoted as R0 Sometimes also called basic reproductive rate or ratio, or basic reproduction ratio Microparasitic infections: R0 is the mean number of secondary infections produced when one infected individual is introduced into a host population where everyone is susceptible Macroparasitic infections: R0 is the mean number of offspring (or female offspring in case of dioecious species) produced by a single parasite, which attain reproductive maturity in a fully susceptible host population

Basic reproduction number

Infectious disease

Czech name

Host

Estimated R0

Measles Pertussis (whooping cough) Chickenpox (varicella) Rubella Mumps Poliomyelitis (polio) Smallpox FIV Rabies Phocine distemper Tuberculosis Influenza Foot-and-mouth disease HIV

Spalničky Černý kašel Plané neštovice Zarděnky Příušnice Dětská obrna Pravé neštovice

Humans (UK) Humans (UK) Humans (UK) Humans (UK) Humans Humans Humans Domestic cats Dogs (Kenya) Seals Cattle Humans Livestock farms (UK) Male homosexuals in England and Wales Female prostitutes in Kenya Humans

16-18 16-18 10-12 6-7 12 5 3.5-6 1.1-1.5 2.44 2-3 2.6 3-4 3.5-4.5 4

HIV Malaria

Vzteklina Chřipka Kulhavka a slintavka

Source: Keeling and Rohani (2008)

11 ≈ 100

Basic reproduction number

R0 is a threshold quantity For most epidemiological models, initial infection Invades a fully susceptible host population if R0 > 1 Dies out if R0 < 1 R0 thus often determines whether or not will an infection invade a host population Exceptions: Seasonally forced models Stochastic models

Basic reproduction number

Simple models with single infected class R0 = product of the Contact rate Mean duration of infection Proportion of infected hosts surviving the latent period of the disease (if there is any)

Basic reproduction number More complex models with several infected classes (e.g. models with age structure) The above simple heuristic definition of R0 for microparasitic infections is insufficient van den Driessche and Watmough (2002) give a precise definition of R0 for a general compartmental epidemiological model defined by ordinary differential equations, and a method of its calculation Definition: R0 is the mean number of secondary infections produced by a typical infected individual in a population at its disease-free equilibrium (DFE) Method of calculation: will be presented later on

Basic reproduction number

Implications for disease control R0 allows us to predict whether an infection invades a host population (R0 > 1) or dies out (R0 < 1) Anything that increases R0 will make invasion and subsequent disease spread more likely Any management strategy (such as vaccination) that decreases R0 below 1 prevents an infection from successfully invading and spreading

Basic reproduction number Replacement number or (effective) reproduction number R For microparasitic infections, it is the mean number of secondary infections produced by a typical infected individual in a host population, not necessarily where all hosts are susceptible or where the population is at its DFE R0 is defined only at the moment of disease invasion, R is defined at all times R0 = R at the moment of disease invasion, as the entire population is susceptible After disease invades, everyone is no longer susceptible: R < R0

Outline of the course

Basic epidemiological models Models of wildlife diseases More complex microparasitic models Seasonality Stochasticity Age structure Spatial dynamics Infectious diseases and control Evolution in host-parasite systems Models of macroparasitic infections

Basic epidemiological models

Agenda for this part Models of microparasitic infections Basic SIR models, incl. models of disease transmission Vaccination Epidemic SIR model with disease-induced mortality Latent period of the disease Infections with a carrier state

Models of microparasitic infections Conventional assumption: host population is held constant, independent of the presence or absence of infection, by an unspecified mechanism Appropriate for many short-term epidemics of non-lethal infections Also for endemic, non-lethal infections of humans in developed countries where natural deaths are roughly balanced by births Mortality induced by lethal infections causes the host population to decrease For endemic diseases, sizes or densities of human populations in developing countries and most animal populations need to be treated as a dynamic variable

Models of microparasitic infections Modelers distinguish several classes of hosts according to their status with respect to the disease

Passively immune (M) IgG antibodies of an infected mother are transferred across the placenta, so that her newborn infants have temporary passive immunity to an infection. M class contains just the infants with passive immunity. After the maternal antibodies disappear from their body, infants move to the susceptible class

Models of microparasitic infections

Modelers distinguish several classes of hosts according to their status with respect to the disease

Susceptible (S) Individuals susceptible to infection. They can contract the disease if they come into contact with it. This class includes also newborn infants who do not have passive immunity, because their mothers were never infected

Models of microparasitic infections

Modelers distinguish several classes of hosts according to their status with respect to the disease

Exposed (E ) Individuals infected but not yet infectious and hence not yet able to pass the disease to others. The period during which individuals are exposed is often called the latent period of the disease

Models of microparasitic infections

Modelers distinguish several classes of hosts according to their status with respect to the disease

Infectious or infective (I ) Individuals that are infectious and hence capable of transmitting the infection to any susceptible that they come into contact with. Do not commute with infected individuals which are those who are either exposed or infectious, E + I

Models of microparasitic infections

Modelers distinguish several classes of hosts according to their status with respect to the disease

Recovered or removed (R) Individuals previously infected but now neither infected nor susceptible. They have an infection-acquired immunity (permanent or temporary)

Models of microparasitic infections Not all epidemiological models include all of these classes, but some include more Choice of which classes to include depends on the characteristics of the modeled disease and the purpose of the model Passively immune class M and exposed class E are often omitted as not considered crucial for the host-parasite interaction Structured models of heterogeneous populations (sex, space, genetic susceptibility, contact pattern, . . .) Acronyms are commonly used to name epidemiological models, based on the classes they contain and the flow patterns between those classes: MSEIR, SEIRS, SIS, . . .

Models of microparasitic infections

Models of microparasitic infections Latent and infectious periods cannot be mistaken with incubation (or asymptomatic) and symptomatic periods, respectively One can transmit parasites long before becoming symptomatic, and one can still be symptomatic while no longer infectious

Source: Anderson and May (1991)

Models of microparasitic infections

Infectious disease

Czech name

Incubation period (day)

Latent period (day)

Infectious period (day)

Measles Mumps Pertussis (whooping cough) Rubella Diphtheria Chicken pox Hepatitis B Poliomyelitis Influenza Smallpox Scarlet fever

Spalničky Příušnice Černý kašel

8-13 12-26 6-10

6-9 12-18 21-23

6-7 4-8 7-10

Zarděnky Záškrt Plané neštovice Žloutenka B Dětská obrna Chřipka Pravé neštovice Spála

14-21 2-5 13-17 30-80 7-12 1-3 10-15 2-3

7-14 14-21 8-12 13-17 1-3 1-3 8-11 1-2

11-12 2-5 10-11 19-22 14-20 1-3 2-3 14-21

Source: Anderson and May (1991)

Models of microparasitic infections

Standard modeling framework Each class = one ordinary differential equation Rate of change of the number or density of individuals in class X : dX = rates of all processes affecting X dt

SIR models Neglect passively immune class M and exposed class E Provide an intuitive basis for understanding disease dynamics dS/dt dI /dt

= =

dR/dt

=

rates of . . . births – natural deaths – new infections rates of . . . new infections – natural deaths – disease-induced deaths – recovery rates of . . . recovery – natural deaths

Model features: Hosts are born as susceptibles S declines due to new infections, I increases at the same rate There are two components of mortality All infectives equally susceptible to disease-induced mortality Immunity is permanent

Epidemic SIR model Epidemic SIR model Good description of rapid disease outbreaks The big question: Why many diseases suddenly develop in a population and then disappear just as suddenly without infecting the entire population? One of the early triumphs of mathematical epidemiology was development of a simple model that was able to predict just this type of behavior observed in countless epidemics Kermack and McKendrick (1927, 1932, 1933)

Epidemic SIR model

How severe will an epidemic be? When will an epidemic start? How many individuals will be affected and thus require treatment? How many individuals will need care at the peak of an epidemic? How long will an epidemic last? How much good would vaccination or another management strategy do in reducing severity of an epidemic? Questions like these can be studied with mathematical models

Epidemic SIR model Kermack-McKendrick epidemic SIR model (core part of almost any epidemiological model) Assumptions: short epidemic ∼ births = natural deaths = 0

non-lethal infection ∼ disease-induced deaths = 0 dS = – rate of new infections dt dI = rate of new infections – rate of recovery dt dR = rate of recovery dt

Modeling disease transmission Description of disease transmission, determining the rate at which new infections occur, is the central part of any epidemiological model and is of crucial importance in determining disease dynamics Key issue – how to model the number of contacts with other individuals per unit time per individual Generally, different types of contact are possible (social, at school, within family, on public transportation, etc.), each occurring at different rate and possibly with different chance per contact of disease being transmitted Kermack-McKendrick epidemic SIR model: homogeneous population = just one type of contacts

Modeling disease transmission

Standard assumption: individuals mix homogeneously (there is no spatial structure) and hence contacts are assumed random Disease transmission rate as a result of random contacts between I and S is a product of four elements: Per individual contact rate with others Φ(N) Proportion of contacts that are with susceptibles S/N Proportion of such contacts that actually result in infection p Number of infectives I rate of new infections = Φ(N)p SI /N

Modeling disease transmission

Rate of occurrence of new cases = number of new cases per unit time = disease incidence Conveys information about the risk of contracting the disease Two most common models of disease transmission: mass action incidence and standard incidence

Modeling disease transmission

Mass action incidence Rate of contacts per individual proportional to population size: Φ(N) = φN Example: doubling the number of people in a bus, a kid with the runny nose is likely to infect twice as many people Φ(N)p SI /N = βSI where β ≡ φp Other names: density-dependent transmission, mass action

Modeling disease transmission Standard incidence Rate of contacts per individual independent of population size: Φ(N) = φ Example: constant number of sexual partners per year or constant number of bites per day by a mosquito Φ(N)p SI /N = βSI /N where β ≡ φp Other names: frequency-dependent transmission, proportionate mixing

Modeling disease transmission

Mass action incidence × standard incidence Much controversy persists over which of these two models is more appropriate (if any) and under what circumstances Mass action incidence (most widely used) Models of air-borne diseases Models of wildlife diseases Standard incidence Models of sexually transmitted diseases Vector-host models

Modeling disease transmission

Mass action incidence × standard incidence Mena-Lorca and Hethcote (1992) rate of new infections = βN v SI /N Data for five human diseases in community sizes between 1,000 and 400,000 leads to values of v between 0.03 and 0.07 Standard incidence (v = 0) likely much more realistic for human diseases than mass action incidence (v = 1) For (these) human diseases the contact rate appears to be largely independent of (or only very weakly dependent on) community size or density

Modeling disease transmission Standard incidence not quite realistic for small population sizes More realistic is to assume that the number of contacts per unit time (i.e. contact rate) per individual is nearly linear for small population sizes (∼ mass action incidence) and saturates for larger population sizes (∼ standard incidence): asymptotic incidence 0.4

0.3 Contact rate

Rate of contacts per individual: N Φ(N) = φ c+N

0.35

0.25 0.2 0.15 0.1 0.05

Φ(N)p SI /N = βSI /(c + N) where β ≡ φp

0 0

1

2 Population density

3

Parameter values: φ = c = 1/3

4

Modeling disease transmission

The assumption of homogeneous mixing: proportion S/N of all contacts made by any I are with S Not always true: e.g. spatial patchiness of infection → variety of other models of disease transmission have been proposed Desperate need for relevant experimental and observational data on transmission dynamics – models of disease transmission and disease dynamics generally outnumber sets of actual data β – transmission coefficient – by far the most difficult to estimate in practical applications

Modeling recovery

Standard assumption: recovery rate scales linearly with the number or density of infectives: rate of recovery = γI γ is the (per capita) recovery rate What does this mean from the perspective of an infectious individual? What is the mean length of the infectious period? Consider a cohort of hosts who are all infectious at time 0, prevent further infections, and let I (τ ) denote the number of hosts who are still infectious at time τ > 0

Modeling recovery With the linear recovery rate, the number of infectives changes as dI /dτ = −γI Solution I (τ ) = I (0)e −γτ ⇒ e −γτ is the proportion of infectives remaining infective at time τ The probability that an individual recovers before time τ after contracting the disease thus is P(recovery time < τ ) = 1 − e −γτ The mean length of the infectious period of an individual: Z ∞ 1 τ γe −γτ dτ = γ 0

Modeling recovery

The linear recovery rate γI thus means that the mean time an individual is infectious is 1/γ and the individual times of infection are exponentially distributed around this mean Not much realistic but mathematically “elegant” More realistic models require a consideration of time delays or a subdivision of the infectious period into smaller intervals and thus become mathematically more complex Many more realistic models exhibit qualitative behavior very similar to the model with linear recovery rate

Epidemic SIR model

Complete epidemic SIR model SI dS = −β(N) dt N dI SI = β(N) − γI dt N dR = γI dt where β(N) ≡ Φ(N)p Initial conditions: a fully susceptible population to which an infective is introduced: S(0) ≈ N, I (0) ≈ 0, R(0) = 0

Epidemic SIR model

Three basic questions regarding epidemics When does an epidemic start? Once an epidemic has started, how will it proceed? How many susceptible individuals will escape an epidemic?

Epidemic SIR model When does an epidemic start? Rate at which new cases are produced by an infectious individual when the entire population is susceptible: β(N)S(0)/N = β(N) Mean infectious period = 1/γ Hence, basic reproduction number R0 = β(N)/γ Infection will spread if initially dI /dt > 0:   S(0) dI = I (0) β(N) −γ >0 dt t=0 N

This happens in a fully susceptible population if β(N)/γ > 1

Hence, an epidemic starts if R0 > 1 and dies out if R0 < 1

Epidemic SIR model When does an epidemic start? Standard incidence: β(N) = β R0 = β/γ Whether an epidemic starts or not depends just on the parasite and host properties Mass action incidence: β(N) = βN R0 = βN/γ It is also the host population size which drives the fate of infection: N > NT ≡ γ/β for an epidemic to start, that is, small enough populations cannot be invaded

Epidemic SIR model Rescaling: divide all three state variables by (constant) host population size N ds = −β(N)si dt di = β(N)si − γi dt dr = γi dt s = S/N, i = I /N, r = R/N are the proportions of hosts in the respective classes Initial conditions: s(0) ≈ 1, i (0) ≈ 0, r (0) = 0

Epidemic SIR model

ds/dt < 0 ⇒ proportion of susceptibles will decrease in time dr /dt > 0 ⇒ proportion of recovered individuals will increase in time i (t) → 0 as t → ∞ Hence, if R0 > 1 then the proportion of infectives will initially increase and eventually decrease to zero Disease spreads, reaches a maximum prevalence and then recedes and vanishes

Epidemic SIR model Once an epidemic has started, how will it proceed?

Population size or density

1000 800 susceptible infectious recovered

600 400 200 0 0

5

10 Time (days)

15

20

Epidemic SIR model Once an epidemic has started, how will it proceed?

Influenza in an English school in 1978, Keeling and Rohani (2008) (N = 763, β = 1.66 per day, 1/γ = 2.2 days: R0 = 3.65)

Epidemic SIR model

Replacement number R = β(N)s(t)/γ decreases as s(t) decreases di = γi (R − 1) dt i will increase for R > 1 or s(t) >

γ β(N)

= 1/R0

i will decline for R < 1 or s(t) < 1/R0 At R = 1 or s = 1/R0 = γ/β(N) disease prevalence reaches a maximum

Epidemic SIR model

Chain rule: di /dt = (di /ds)(ds/dt) ⇒  ds di γ di = = −1 + ds dt dt β(N)s i = −s +

γ γ ln s + c, c = i (0) + s(0) − ln s(0) β(N) β(N)

Inserting s = 1/R0 , s(0) = 1 and i (0) = 0 implies imax = 1 −

1 1 1 + ln R0 R0 R0

Epidemic SIR model Once an epidemic has started, how will it proceed? 1000

β = 1 (R0 = 2) β = 2.5 (R0 = 5)

Disease incidence

800

β = 4 (R0 = 8) 600 400 200 0 0

5

10 Time

15

20

Higher values of R0 lead to a shorter yet more severe epidemic, as dimax /dR0 > 0 (N = 1000, γ = 0.5)

Epidemic SIR model How many susceptibles will escape an epidemic? Chain rule, ds/dt = (ds/dr )(dr /dt) ⇒  ds dr β(N) ds = =− s = −R0 s dr dt dt γ r (0) = 0 ⇒ s(t) = s(0) exp(−r (t)R0 ) As r (t) ≤ 1 ⇒ s(t) ≥ s(0) exp(−R0 ) ⇒ there will always be susceptibles that will escape epidemic s + i + r = 1 and i (t) → 0 as t → ∞ ⇒ s∞ = s(0) exp(−(1 − s∞ )R0 ) where s(t) → s∞ as t → ∞ is the final proportion of susceptibles

Epidemic SIR model How many susceptibles will escape an epidemic? 1

s∞ (blue), r∞ (red)

0.8 0.6 0.4 0.2 0 1

2

3

4 R0

5

Higher values of R0 lead to a lower s∞

6

7

Estimating R0 for an epidemic

Let an isolated village experience an outbreak of influenza in which 812 of the 1100 residents contract the infection. Estimate R0 assuming that the outbreak started with a single case contracted from outside the willage, with all others susceptible at the start of the outbreak. s∞ = s(0) exp(−(1 − s∞ )R0 ) implies R0 = − We know that r∞ = 812/1100 and hence s∞ = 1 − r∞ = 288/1100, and that s(0) = 1 As a result, R0 ≈ 1.8

ln(s∞ /s(0)) 1 − s∞

Endemic SIR model

Because of permanent immunity to reinfection, disease epidemics race through the host population as waves of infection New susceptible individuals are born behind these waves To model an endemic disease one needs to think of a longer time scale and thus include births and deaths The simplest way: birth rate = natural death rate, no deaths due to the disease ⇒ total population size remains constant

Endemic SIR model Model equations SI dS = µN − β(N) − µS dt N SI dI = β(N) − γI − µI dt N dR = γI − µR dt All hosts produce offspring at the same rate (equal to the natural mortality rate) µ; 1/µ can thus be viewed as the natural host life expectancy Both birth and death rates are assumed density-independent

Endemic SIR model Threshold behavior Basic reproduction number R0 =

β(N) γ+µ

Invasion criterion (dI /dt)|t=0,S=N = I (β(N) − γ − µ) > 0 Disease invades if R0 > 1 and dies out if R0 < 1 R0 lower here than for the epidemic SIR model – death rate reduces the mean length of the period in which an individual is infectious

Endemic SIR model Rescaling: divide all three state variables by N: ds/dt = µ − β(N)si − µs

di /dt = β(N)si − γi − µi dr /dt = γi − µr

Disease-free equilibrium (DFE) (1, 0, 0) Endemic equilibrium (s ∗ , i ∗ , r ∗ ) with i ∗ > 0: 1 µ s∗ = (R0 − 1), r ∗ = 1 − s ∗ − i ∗ , i∗ = R0 β(N) Endemic equilibrium is feasible (0 < i ∗ < 1) if and only if R0 > 1, which makes perfect sense as we know that the disease can invade and persist in the host population iff R0 > 1

Endemic SIR model Local asymptotic stability of model equilibria Jacobian: J(s, i ) =

"

−β(N)i − µ −β(N)s

−β(N)i

β(N)s − γ − µ

#

Disease-free equilibrium (1, 0): " # −µ −β(N) J(1, 0) = 0 β(N) − γ − µ Locally asymptotically stable ⇔ β(N) − γ − µ < 0 ⇔ R0 < 1, otherwise unstable

Endemic SIR model

Local asymptotic stability of model equilibria Endemic equilibrium (s ∗ , i ∗ ): " # −µR0 −(γ + µ) ∗ ∗ J(s , i ) = µ(R0 − 1) 0 TrJ = −µR0 < 0

det J = µ(R0 − 1)(γ + µ) > 0 ⇔ R0 > 1 Locally asymptotically stable ⇔ R0 > 1 ⇔ it exists

Endemic SIR model

Local stability analysis If R0 < 1 (disease dies out), the DFE is locally asymptotically stable If R0 > 1 (disease invades) the DFE is unstable and the endemic equilibrium is locally asymptotically stable R0 > 1: supplementing the pool of susceptibles by newborns ensures disease persistence in the long run

Endemic SIR model Solutions might approach the endemic equilibrium in an oscillatory manner, with the amplitude declining over time

Size or density of infectives

2

1.5

1

0.5

0 0

10

20

30 40 Time (years)

50

60

Endemic SIR model

Mean age at infection Force of infection: The per capita rate at which susceptibles catch the disease = the probability per unit time that a susceptible individual becomes infected Force of infection at endemic equilibrium: F = βI ∗ /N Ignoring small disease-independent mortality, the inverse of F is an approximate mean time an individual spends in the susceptible class = mean age at infection A if all individuals are born susceptible

Endemic SIR model

Mean age at infection Standard incidence: F = βI ∗ /N and I ∗ = N(µ/β)(R0 − 1): A = 1/F ≈

1 µ(R0 − 1)

This can be rephrased as R0 ≈ 1 + L/A where L = 1/µ is the host life expectancy

Endemic SIR model

Mean age at infection (years)

80 1/µ = 40 years 1/µ = 60 years 1/µ = 80 years

60

40

20

0 2

4

6

8

10

12 R0

14

16

18

20

Endemic SIR model

Mean age at infection (in years) for different diseases in different places and periods (Anderson and May 1991)

SIS model Endemicity is strengthened when there is no or only temporary immunity – susceptible class is then steadily supplied by individuals that recover from infection SIS model dS SI = µN − β(N) − µS + γI dt N SI dI = β(N) − µI − γI dt N SIS models are commonly used to study sexually transmitted diseases, where due to the vast number of strains there is no such thing as a recovered class

SIS model

1000 Size or density infectious

susceptible infectious 800 600 400 200 0 0

20

40 60 Time (days)

80

100

Temporal disease dynamics when the disease attains the endemic equilibrium (N = 1000; µ = 0.1 per day, β = 1.66 per day, 1/γ = 2.2 days)

SIRS model Endemicity is strengthened when there is no or only temporary immunity – susceptible class is then steadily supplied by individuals that recover from infection SIRS model SI dS = µN − β(N) − µS + σR dt N SI dI = β(N) − γI − µI dt N dR = γI − µR − σR dt σ = rate at which immunity is lost and recovered individuals return to the susceptible class; 1/σ can thus be viewed as the mean duration of immunity against reinfection

SIRS model

1 susceptible infectious recovered

800

Equilibrial infection prevalence i*

Population size or density

1000

600 400 200 0 0

10

20 30 Time (days)

40

50

0.8 0.6 0.4 0.2 0 0

1 2 3 4 Length of immunity period 1/σ

5

Left: Temporal disease dynamics when the disease attains the endemic equilibrium (N = 1000; µ = 0.1 per day, β = 1.66 per day, 1/γ = 2.2 days, 1/σ = 3.2 days) Right: The effect of waning immunity on the equilibrium infection prevalence (µ = 0.1 per day, β = 1.66 per day, 1/γ = 2.2 days)

Summary

Host population size/desity constant; no disease-induced mortality Epidemic models: disease extinction or a wave-like epidemic followed by disease extinction Endemic models: disease extinction or an approach (monotonous or oscillatory) to an endemic equilibrium For a disease to be endemic, there must be a supply of fresh susceptibles: through births (epidemic SIR model) plus through recovery with no immunity (SIS model) or temporary immunity (SIRS model)

SIR models

Limitations to simple SIR models Host population size is constant, i.e. no disease-induced mortality is considered There is no latent period of the disease There are no dynamics of the host population when the disease is absent The host population is uniform and homogeneously mixing (in reality, children usually have more contacts per day than adults, different geographic and socio-economic groups have different contact rates, etc.)

SIR models

Main value of simple SIR models is in letting us know what is possible in an as simple and abstract system as possible where everything except the actual host-parasite interaction has been removed So, despite these limitations, these simple SIR models have contributed much to our understanding of the course of real epidemics and the way endemic equilibria can be attained, and capture much of the qualitative properties of more complex, detail-rich models In what follows, we will see models which are free of some of these limitations

Vaccination

Herd immunity = enough people have disease-acquired or vaccination-acquired immunity, so that introduction of one infective into the population does not cause an invasion of the disease The ultimate goal of any vaccination program = vaccinate enough people to achieve herd immunity Successful vaccination moves individuals from the susceptible class straight into the recovered class so that they can no longer catch or spread the infection

Vaccination Vaccination of newborns Let p be the proportion of newborns vaccinated: dS SI = µN(1 − p) − β(N) − µS dt N SI dI = β(N) − γI − µI dt N dR = µNp + γI − µR dt Disease cannot invade the host population if p > pc = 1 − (γ + µ)/β(N) = 1 − 1/R0

Vaccination Disease

R0

pc [%]

Smallpox Measles Chickenpox

3-5 16-18 8-10

66-80 93-95 87.5-90

Source: Keeling and Rohani (2008) Smallpox: the only infectious disease for which successful vaccination has been achieved worldwide and smallpox has been eliminated (the last known case was in Somalia in 1977) Measles: herd immunity against measles has not been achieved and probably will never be, as the minimum required proportion of newborns vaccinated is extremely high (and the vaccine for measles is not always effective)

Vaccination

Vaccination of newborns If p < pc and the disease persists, the equilibrium disease prevalence (i.e. proportion of infectious individuals) is: i∗ =

I∗ µ(1 − p) µ = − N γ+µ β(N)

Thus, vaccination reduces the disease prevalence

Vaccination Vaccination of susceptibles Let θ be the rate at which susceptibles are vaccinated (vaccination of I and R individuals is assumed to have no effect): SI dS = µN − β(N) − µS − θS dt N SI dI = β(N) − γI − µI dt N dR = θS + γI − µR dt Disease cannot invade the host population if θ > θc = µ(β(N)/(γ + µ) − 1) = µ(R0 − 1)

Vaccination

In reality, vaccination is only partly effective, decreasing the rate of new infections and also decreasing infectivity of a vaccinated person which has become infected. More complex models are needed to come up with sound quantitative recommendations on applicable vaccination strategies.

Why do we care so much about R0 anyway? R0 provides five fundamental insights into the dynamics of an infectious disease: R0 is the threshold parameter, determining whether or not there will be an epidemic R0 determines the initial rate of increase of an epidemic (i.e., during its exponential growth phase) R0 determines the final size of the epidemic (i.e., what fraction of susceptibles will ultimately be infected over the course of the outbreak) R0 determines the endemic equilibrium proportion of susceptibles in the population R0 determines the critical vaccination threshold

Epidemic SIR model with disease-induced mortality Epidemic SIR model Assumptions: Infectious individuals leave the I class at rate αI Proportion f of these individuals recover Remaining proportion 1 − f of them die of the disease dS/dt = −Φ(N)p SI /N

dI /dt = Φ(N)p SI /N − αI dR/dt = f αI

Total host population size is no more constant: dN/dt = −(1 − f )αI

Epidemic SIR model with disease-induced mortality

All hosts are initially susceptible, N(0) = K ⇒ R0 = Φ(K )p/α Since I 0 (0) = α(R0 − 1)I (0), then if R0 > 1 an epidemic starts and if R0 < 1 the disease dies out

0

0

S + I = −αI ⇒ S(t) + I (t) − S(0) − I (0) = −α t→∞⇒α

Z

0



I (s)ds = K − S∞

Z

0

t

I (s)ds

Epidemic SIR model with disease-induced mortality

dS/dt = −Φ(N)p SI /N ⇒ −

S0 Φ(N) = pI S N

Let Φ(N) > 0, Φ0 (N) > 0, and (Φ(N)/N)0 ≤ 0 (saturating) Let B = limN→0 Φ(N)/N be finite ⇒ Φ(N)/N ≤ B Integration from 0 to ∞ gives Z ∞ Z ∞ Φ(N(t)) Bp(K − S∞ ) S(0) I (t)dt = pI (t)dt ≤ Bp = ln S∞ N(t) α 0 0 Since the right-hand side is finite, the left-hand side is also finite and hence S∞ > 0

Epidemic SIR model with disease-induced mortality The epidemic SIR model with disease-induced mortality demonstrates the same qualitative behavior as the simple epidemic SIR model with constant population size: R0 = Φ(K ) p/α decides on disease invasion Some hosts escape the disease after the epidemic passes 1000 Population size or density

Final susceptible Final recovered 800 600 400 200 0 0

0.2

0.4 0.6 Fraction recovering f

0.8

1

Parameter values: Φ(N) = bN/(1 + cN), b = 1, c = 3, p = 4, α = 0.5, N(0) = K = 1000

Latent period: SEIR model

Infectious disease

Czech name

Incubation period (day)

Latent period (day)

Infectious period (day)

Measles Mumps Pertussis (whooping cough) Rubella Diphtheria Chicken pox Hepatitis B Poliomyelitis Influenza Smallpox Scarlet fever

Spalničky Příušnice Černý kašel Zarděnky Záškrt Plané neštovice Žloutenka B Dětská obrna Chřipka Pravé neštovice Spála

8-13 12-26 6-10 14-21 2-5 13-17 30-80 7-12 1-3 10-15 2-3

6-9 12-18 21-23 7-14 14-21 8-12 13-17 1-3 1-3 8-11 1-2

6-7 4-8 7-10 11-12 2-5 10-11 19-22 14-20 1-3 2-3 14-21

Source: Anderson and May (1991)

Latent period: SEIR model Exposed class E : parasites reproduce rapidly within a host but their abundance is still too low for a transmission to other susceptible hosts dS/dt = µN − β(N) SI /N − µS dE /dt = β(N) SI /N − σE − µE dI /dt = σE − γI − µI dR/dt = γI − µR σ = rate at which individuals leave the exposed class and enter the infectious one; 1/σ can thus be viewed as the mean latent period σ → ∞: latent period negligible, SEIR model → SIR model

Latent period: SEIR model Basic reproduction number R0 =

β(N) σ × γ+µ µ+σ

is the product of the contact rate β(N) mean proportion surviving the latent period σ/(µ + σ) mean length of the infectious period 1/(γ + µ) An individual has to survive the latent period in order to produce new infectives Typically σ/(µ + σ) ∼ 1 (latent period usually  mean life expectancy, e.g. humans in developed countries have 1/µ ∼ 70 years): SEIR model → SIR model

Latent period: SEIR model Dynamic properties of the SEIR model are qualitatively analogous to those of the endemic SIR model R0 < 1 ⇒ DFE locally asymptotically stable R0 > 1 ⇒ DFE unstable and a unique endemic equilibrium exists and is locally asymptotically stable

Equilibrium size or density

3

x 10

B

−3

susceptible exposed infectious

2.5

Equilibrium size or density

A

2 1.5 1 0.5 0 0

2

4 6 8 Latent period 1/σ (days)

10

0.5 susceptible exposed infectious

0.4 0.3 0.2 0.1 0 0

2

4 6 8 Latent period 1/σ (days)

10

Parameter values: N = 1000, β = 73 per day, 1/γ = 7 days, 1/µ = 70 years (A), 1/µ = 7 days (B)

Latent period: SEIR model

dS/dt = −β(N) SI /N dE /dt = β(N) SI /N − σE dI /dt = σE − γI dR/dt = γI

Number of infectious individuals

400

Epidemic SEIR model

1/σ = 0 1/σ = 1 days 1/σ = 2 days

350 300 250 200 150 100 50 0 0

10

20 Time (days)

30

40

Parameter values: N = 1000; β = 1.66, γ = 1/2.2

Longer latent periods cause the epidemic to begin at a slower rate, reach a lower peak prevalence, but last much longer

Infections with a carrier state: SICR model Hepatitis B or herpes: a proportion of infected individuals may become chronic carriers, transmitting infection at a low rate for many years after recovering We focus on development of a model with a single carrier class Susceptible individuals can be infected by either carriers of a novel class C or by infectious individuals of the infectious class I Progress of infection within an individual is assumed independent of the source of infection (whether C or I ) Recently infected individual is infectious for a given period and then either recovers completely or moves into the carrier class

Infections with a carrier state: SICR model Model

dS/dt = µN − β SI /N − β SC /N − µS dI /dt = β SI /N + β SC /N − γI − µI dC /dt = qγI − ΓC − µC dR/dt = (1 − q)γI + ΓC − µR

C = abundance of the carriers in the population  = reduced transmission rate from chronic carriers compared to infectious individuals (0 <  < 1) q = proportion of infectious individuals that become carriers while a fraction 1 − q simply recover (0 < q < 1) Γ = rate at which individuals leave the carrier class to eventually recover

Infections with a carrier state: SICR model

Basic reproduction number R0 = sum of two components, one coming from the infectious individuals and the other from the chronic carriers R0 =

qγ β β + × γ+µ γ+µ Γ+µ

qγ/(γ + µ) = fraction of those individuals that do not die when in the infectious class and go on to become chronic carriers The fact that infectious individuals can enter the carrier state rather than simply recover increases the value of R0 and hence the likelihood of disease invasion

Infections with a carrier state: SICR model

A

B 0.05

0.06 susceptible infectious carriers

Proportion of population

Proportion of population

0.06

0.04 0.03 0.02

susceptible infectious carriers

0.04 0.03 0.02 0.01

0.01 0 0

0.05

0.2 0.4 0.6 0.8 Fraction becoming carrier, q

1

0 0

0.2 0.4 0.6 0.8 Carriers transmission reduction factor, ε

1

Parameter values: β = 73 per year, 1/γ = 100 days, 1/Γ = 1000 days, 1/µ = 50 years,  = 0.1 (A), q = 0.3 (B)

Modeling a zombie infection

Munz et al. (2009) When zombies attack!: Mathematical modelling of an outbreak of zombie infection. In: Infectious Disease Modelling Research Progress, pp. 133-150. Nova Science Publishers, Inc. Susceptibles (S), zombie (Z ), removed (R) dS/dt = Π − β SZ − δS dZ /dt = β SZ + ηR − α SZ dR/dt = δS + α SZ − ηR

Modeling a zombie infection

dS/dt + dZ /dt + dR/dt = Π S + Z + R → ∞ as t → ∞, if Π 6= 0 Equilibrium analysis: S ∗ =

Π ⇒S 9∞ βZ ∗ + δ

Doomsday scenario: an outbreak of zombies will lead to the collapse of civilization, as large numbers of people are either zombified or dead

Modeling a zombie infection

If an outbreak occurs on a short timescale, we can assume Π=δ=0 N = population size Disease-free equilibrium: (N, 0, 0) unstable Doomsday equilibrium: (0, Z , 0) locally stable Adding a latent period: the same qualitative outcome

Modeling a zombie infection

A model with treatment dS/dt = Π − β SZ − δS + cZ dI /dt = β SZ − ρI − δI dZ /dt = ρI + ηR − α SZ − cZ dR/dt = δS + δI + α SZ − ηR Π 6= 0 = doomsday scenario, we again assume Π = δ = 0

Modeling a zombie infection

Disease-free equilibrium: (N, 0, 0, 0) unstable Coexistence equilibrium:   ∗ locally stable (S ∗ , I ∗ , Z ∗ , R ∗ ) = βc , ρc Z ∗ , Z ∗ , αc Z ηβ Humans are not eradicated, but only exist in low numbers

Modeling wildlife diseases

Agenda for this part Why to study wildlife diseases? SI models with exponential host population growth SI models with logistic host population growth SI models with constant rate immigration of hosts

Why to study wildlife diseases? Three major reasons for studying wildlife diseases Wildlife species of conservation interest – concerns about the impact of diseases on population survival spillover from domestic animals: canine distemper virus from domesticated dogs to African wild dogs diseases to only impact wildlife: chytridiomycosis (fungal infection) in amphibians

Increasing concerns about the possibility of either transmission from wildlife to humans or to domestic animal species – wildlife as the reservoir species: HIV, Ebola, bird flu, . . . Many natural phenomena driven by parasites (sexual selection, population dynamics, . . .)

Why to study wildlife diseases?

Science, Vol. 287, No. 5452 (Jan. 21, 2000), pp. 443-449

Why to study wildlife diseases? Differences between modeling human and wildlife diseases Wildlife populations do not remain constant over time. Populations can be highly variable due to environmental factors, landscape, or because of their internal dynamics. Multiple species interactions are often involved. In many cases, wildlife populations are believed to be controlled by parasites. Parasites can thus be added to the weaponry aimed at controlling pest populations, or limited in their ability to invade endangered populations of hosts, to prevent biodiversity degradation. Data can be more easily obtained from animal disease systems, since it is often possible to experiment on wildlife populations or individual animals without the ethical issues involved with human disease systems.

SI model Generic SI model that we will study dS/dt = b(N)N − β(N)SI /N − d(N)S

dI /dt = β(N)SI /N − d(N)I − αI

b(N) and d(N) = possibly density-dependent birth and death rates Total host population N = S + I evolves as dN/dt = b(N)N − d(N)N − αI Mass action incidence: β(N) = βN Standard incidence: β(N) = β

Exponential host population growth Host population in the absence of disease is assumed to grow exponentially Rapidly multiplying animals such as insects (especially at invasion or outbreak) 20

dN/dt = bN − dN

Solution: N(t) = N(0)e rt r =b−d

Population size or density

b(N) = b, d(N) = d

r=1 r=2 r = −1

15

10

5

r > 0 . . . exponential growth r < 0 . . . decline to extinction

0 0

0.5

1 Time

1.5

2

Exponential host population growth

Major question Whether and how a disease may modify dynamics of an exponentially growing host population Infections are expected to slow down population growth Will a disease regulate the host population to a steady state or even make it extinct?

Exponential host growth & mass action incidence

SI model with mass action incidence dS = bN − βSI − dS dt dI = βSI − dI − αI dt Total host population evolves as dN/dt = rN − αI where r = b − d Assumption: b > d

Exponential host growth & mass action incidence

Basic reproduction number βN(0) βS = R0 = α + d S=N(0) α+d

Disease invades a fully susceptible host population if R0 > 1 R0 > 1 if and only if N(0) > NT ≡ (α + d)/β N(0) > NT ⇒ I will increase → disease will spread N(0) < NT ⇒ I will decline and S will grow at rate ∼ r until N > NT ; I will then increase and disease will spread

Exponential host growth & mass action incidence

Equation for I dI /dt = βI (N − I − NT ) Candidates for equilibria N∗ =

α r NT and I ∗ = NT ⇒ S ∗ = NT α−r α−r

I ∗ > 0 provided that α > αT ≡ r This positive equilibrium is unique and locally stable

Exponential host growth & mass action incidence

Relatively mild disease α < αT No positive equilibrium exists, host population grows exponentially dN = rN − αI = (r − α)N + αS > (r − α)N dt This implies N(t) > N(0)e (r −α)t

Exponential host growth & mass action incidence

Dynamics of i = I /N (recall N is a function of time) di /dt = i [βN(1 − i ) − (α + d + r ) + αi ] N grows exponentially and 0 ≤ i ≤ 1 ⇒ i → 1 dN/dt = N(r − αi ) → N(r − α) ⇒ N grows at rate ρ = r − α N(t) → ce ρt , i → 1 and I = iN ⇒ I (t) → ce ρt dI /dt = βI (S − NT ) ⇒ ρ = β(Sˆ − NT ) ⇒ S → Sˆ = b/β

Exponential host growth & mass action incidence

40 35

no disease disease, α > α

30

disease, α < α

Total host population

T T

25 20 15 10 5 0 0

50

100

150 Time

200

250

300

Exponential host growth & standard incidence

SI model with standard incidence dS SI = bN − β − dS dt N SI dI = β − dI − αI dt N Total host population N = S + I evolves as dN/dt = rN − αI where r = b − d Assumption: b > d

Exponential host growth & standard incidence

Basic reproduction number βS/N(0) β R0 = = α + d S=N(0) α+d

Disease invades a fully susceptible host population if R0 > 1 R0 > 1 . . . disease invades R0 < 1 . . . disease declines to extinction; host population goes on to grow at rate r

Exponential host growth & standard incidence

Dynamics of i = I /N di /dt = d(I /N)/dt = (dI /dt − i (dN/dt))/N di = i [β(1 − i ) − (α + d + r ) + αi ] dt Proportional equilibria: i0∗ = 0 and ie∗ = 1 − b/(β − α) i0∗ always exists ˆ ≡ β/(b + α) > 1 ie∗ is feasible (0 < ie∗ < 1) if R

Exponential host growth & standard incidence

i0∗ = 0 Locally stable if unique, otherwise unstable ˆ < 1), dN/dt = (r − αi )N → rN: host population Once stable (R grows exponentially at rate r dI /dt = I [βS/N − (d + α)] → I (d + α)(R0 − 1) as t → ∞ R0 > 1 . . . disease invades and I → ∞ R0 < 1 . . . disease dies out and I → 0

Exponential host growth & standard incidence

i ∗ = 1 − b/(β − α) ˆ > 1: existence of i ∗ > 0 implies disease invasion b > d ⇒ R0 > R e Locally stable if it exists dN/dt = N(r − αie∗ ): N → 0 if r − αie∗ < 0 or R1 ≡ b/(d + αie∗ ) < 1 N → ∞ exponentially at rate ρ = r − αie∗ if R1 > 1

Exponential host growth & standard incidence

40

Total host population

35 30

ˆ 1, R1 > 1 ˆ > 1, R1 < 1 disease, R

25 20 15 10 5 0 0

20

40

60 Time

80

100

SI model with exponential host population growth Qualitatively different results 40 35

40 no disease disease, α > α

35

disease, α < α

30

30

Total host population

Total host population

T T

25 20 15 10 5 0 0

ˆ 1, R1 > 1 disease, R ˆ > 1, R1 < 1 disease, R

25 20 15 10 5

50

100

150 Time

200

250

300

0 0

20

40

60

80

Time

Mass action incidence:

Standard incidence:

The host population

The host population

grows exp. at reduced rate

grows exp. at reduced rate

attains a stable eq.

goes extinct

100

Logistic host population growth Host population in the absence of disease is assumed to grow logistically Many populations are limited in growth by a finite amount of resources (such as food or territories) – the per capita population growth rate then declines as population density increases 8

r > 0 . . . intrinsic growth rate K > 0 . . . carrying capacity

N(0) = 0.1 N(0) = 3 N(0) = 7

7 Population size or density

Logistic equation:   dN N = rN 1 − dt K

6 5 4 3 2 1 0 0

2

4

6 Time

8

10

Logistic host population growth Drawback: no distinction between birth and death rates Negative density dependence: birth rate decreases and/or death rate increases with increasing population density Birth rate: b(N) = b/(1 + b1 N), b(N) = b exp(−b1 N), or b(N) = b − b1 N if N ≤ b/b1 and b(N) = 0 otherwise Death rate: d(N) = d + d1 N, d(N) = dN/(d1 + N) d . . . intrinsic death rate d1 . . . strength of negative density-dependent effects dN = bN − (d + d1 N)N dt

Logistic host growth & mass action incidence SI model with mass action incidence dS = bN − βSI − (d + d1 N)S dt dI = βSI − (d + d1 N)I − αI dt Equal density dependence for susceptibles and infectives Total host population N = S + I evolves as dN/dt = bN − (d + d1 N)N − αI = rN(1 − N/K ) − αI where r = b − d and K = (b − d)/d1 Assumption: b > d

Logistic host growth & mass action incidence Basic reproduction number βK βS = R0 = d + d1 N + α S−N=K d + d1 K + α =

βK βK = d + d1 ((b − d)/d1 ) + α b+α

Disease invades a fully susceptible host population at carrying capacity (i.e. at DFE) if R0 > 1 ⇔ β > (b + α)/K If R0 < 1 or β < (b + α)/K , the disease will die out ⇒ populations at low K cannot be invaded

Logistic host growth & mass action incidence

Let us analyze equations for i = I /N and N: di = i [β(1 − i )N − α − b + αi ] dt dN = N[b − (d + d1 N) − αi ] dt Four candidate equilibria: (0, 0), (0, K = (b − d)/d1 ), (1 + b/α, 0), and endemic equilibrium (ie∗ , Ne∗ ) with ie∗ > 0 and Ne∗ > 0 Equilibrium (1 + b/α, 0) not feasible as i ∗ = 1 + b/α > 1

Logistic host growth & mass action incidence Endemic equilibrium (ie∗ , Ne∗ ) exists provided that R0 > 1 It solves β(1 − i )N − α − b + αi = 0 and b − (d + d1 N) − αi = 0

Logistic host growth & mass action incidence

Local stability of model equilibria Jacobian:  J =

−2βNi + 2αi + βN − α − b β(1 − i )i −αN

−2d1 N + b − d − αi

Endemic equilibrium (ie∗ , Ne∗ ): TrJ = −i (βN − α) − d1 N det J = d1 Ni (βN − α) + αNβ(1 − i )i

 

Logistic host growth & mass action incidence

Endemic equilibrium (ie∗ , Ne∗ ): Sufficient condition for local stability: βN > α At endemic equilibrium, the replacement number R = 1: R=

βN(1 − i ) βS = =1 d +α d +α

βN − α = d + i βN > 0 βN > α

Logistic host growth & mass action incidence

Total host population

12 10 8 6 4 no disease; disease, R0 < 1

2 0 0

disease, R0 > 1 200

400

600 Time

800

1000

Logistic host growth & standard incidence SI model with standard incidence SI dS = bN − β − (d + d1 N)S dt N SI dI = β − (d + d1 N)I − αI dt N Equal density dependence for susceptibles and infectives Total host population N = S + I evolves as dN/dt = bN − (d + d1 N)N − αI = rN(1 − N/K ) − αI where r = b − d and K = (b − d)/d1 Assumption: b > d

Logistic host growth & standard incidence

Basic reproduction number βS/N β R0 = = d + d1 N + α S=N=K b+α

Disease invades a fully susceptible host population at carrying capacity (i.e. at DFE) if R0 > 1 ⇔ β > b + α If R0 < 1 or β < b + α, the disease will die out Populations of any size can be invaded

Logistic host growth & standard incidence

Dynamics of s = S/N and i = I /N ds = b(1 − s) + (α − β)si dt di = βsi − i (α + b) + αi 2 dt s =1−i di = i [β(1 − i ) − (α + b) + αi ] dt i ∗ = 0 ⇒ s ∗ = 1 ⇒ N → K . . . disease-free equilibrium (DFE)

Logistic host growth & standard incidence Endemic proportional equilibrium: i∗ = 1 −

b β−α

and s ∗ =

b β−α

0 < i ∗ < 1 as soon as R0 > 1 dN/dt = N[b − (d + d1 N) − αi ] Two equilibria for N: N ∗ = 0 and N ∗ = (b − d − αi ∗ )/d1 Local stability analysis: N → N ∗ as t → ∞ where N ∗ = K if R0 < 1

N ∗ = (b − d − αie∗ )/d1 < K if R0 > 1 and b − d > αie∗ N ∗ = 0 if R0 > 1 and b − d < αie∗

Logistic host growth & standard incidence

Total host population

12 10 8 6

no disease; disease, R0 < 1

4

disease, R0 > 1, i* < (b−d)/α

2

disease, R0 > 1, i* > (b−d)/α

0 0

500

1000 Time

1500

2000

SI model with logistic host population growth Almost identical results The host population attains a stable equilibrium which might be reduced due to the disease 12 Total host population

Total host population

12 10 8 6 4 no disease; disease, R0 < 1

2 0 0

disease, R0 > 1 200

400

600 Time

Mass action incidence: Disease cannot make the population extinct

800

1000

10 8 6

no disease; disease, R0 < 1

4

disease, R0 > 1, i* < (b−d)/α

2

disease, R0 > 1, i* > (b−d)/α

0 0

500

1000 Time

1500

Standard incidence: Disease can make the population extinct

2000

Constant rate immigration of hosts Many models assume immigration at a constant rate instead of a density-dependent birth term Immigration rate Λ > 0 8 N(0) = 0.1 N(0) = 3 N(0) = 7

dN = Λ − dN dt

Population size or density

7 6 5 4 3 2 1 0 0

2

4

6

8

10

Time

Host population never goes extinct in the absence of disease and always attains the carrying capacity Λ/d

Constant rate immigration of hosts

SI model with generic incidence SI dS = Λ − β(N) − dS dt N SI dI = β(N) − dI − αI dt N Total host population N = S + I evolves as dN/dt = Λ − dN − αI

Constant rate immigration of hosts Disease-free equilibrium: (Λ/d, 0) For the disease to invade the host population we must have 1 dI = β(N)S/N − d − α > 0 I dt when evaluated at DFE ⇒ β(Λ/d) > d + α Basic reproduction number: β(Λ/d) β(N)S/N = R0 = d + α S=N=K =Λ/d d +α

⇒ disease invades if R0 > 1

DFE locally stable if R0 < 1 and unstable if R0 > 1

Constant rate immigration of hosts

Endemic equilibrium must satisfy: β(N)S = (d + α)N (from the equation for I ) Λ + αS = (d + α)N (from the equation for N) This can be solved for S and N if a specific form of β(N) is given Once solved, the equilibrium density of infectives is I =

Λ − dS (from the equation for S) (d + α)

Constant rate immigration of hosts

Standard incidence β(N) = β: a unique solution S∗ =

β − (d + α) Λ and I ∗ = Λ β−α (d + α)(β − α)

Mass action incidence β(N) = βN: a unique solution S∗ =

d +α Λβ − d(d + α) and I ∗ = β β(d + α)

Both cases: S ∗ > 0 and I ∗ > 0 (feasible) if R0 > 1 This unique endemic equilibrium is locally stable if it exists

SI model with constant rate immigration of hosts Identical results The host population always attains a stable equilibrium which might be reduced due to the disease The host population can never decline to zero Mass action incidence:

Standard incidence: 20 Total host population

Total host population

20

15

10

5

no disease; disease, R0 < 1

15

10

5

no disease; disease, R0 < 1

disease, R0 > 1 0 0

50

100 Time

150

disease, R0 > 1 200

0 0

50

100 Time

150

200

Overall summary

Exponential growth

Logistic growth

Mass action incidence

Mass action incidence

Constant immigration Mass action incidence

no disease disease, α > αT

30

disease, α < α

20

12

T

25 20 15 10

Total host population

35

Total host population

Total host population

40

10 8 6 4 no disease; disease, R0 < 1

2

5

15

10

5

no disease; disease, R0 < 1

disease, R > 1

disease, R > 1

0

50

100

150 Time

200

250

0 0

300

30

600

0

800

0 0

1000

Standard incidence

ˆ 1, R1 > 1 disease, R ˆ > 1, R1 < 1 disease, R

50

25 20 15 10

150

200

Standard incidence

10 8 6

no disease; disease, R0 < 1

4

disease, R > 1, i* < (b−d)/α 0

disease, R > 1, i* > (b−d)/α

2

5

100 Time

20

12 Total host population

Total host population

35

400 Time

Standard incidence 40

200

Total host population

0 0

15

10

5

no disease; disease, R0 < 1

0

disease, R > 1 0

0 0

20

40

60 Time

80

100

0 0

500

1000 Time

1500

2000

0 0

50

100 Time

150

200

Overall summary

Several important points Parasite regulation may involve a qualitative change. Exponentially growing populations can be stabilized or even made extinct, just like populations that grow logistically in the absence of infection. Disease-induced extinction requires the transmission rate not to decline to zero as population density declines (true for standard incidence). Otherwise, low-density populations are virtually free of parasites and can recover.

Overall summary

Several important points The exception is models with constant immigration rate. Since in this case susceptibles arrive at a steady rate, such populations cannot be made extinct by any disease no matter what incidence drives parasite transmission. The basic types of dynamics revealed for our simplest SI models are conserved in more complex epidemiological models.

Overall summary A number of decisions are required when considering a host-parasite interaction, even in the framework of the simplest host-parasite models, including type of host demography (constant population density [epidemic or endemic disease], exponential growth, logistic growth, or constant rate immigration) whether to consider no, temporary or permanent immunity to reinfection which form to use for (horizontal) transmission rate (frequency-dependent, density-dependent or other) if there is no or some disease-induced mortality if there is no or some vertical transmission if there is no or some depression of reproduction due to infection

Overall summary

Simple extension for infectious diseases that reduce host fecundity and/or are transmitted vertically: SI dS = b(S + (1 − p)(1 − σ)I ) − β(N) − dS dt N SI dI = bp(1 − σ)I + β(N) − dI − αI dt N σ . . . proportional decrease in fecundity of infectives p . . . proportion of newborns produced by infectives to which the disease is transmitted vertically

Overall summary

Essential message Qualitative properties at the level of individual hosts and parasites might create qualitative differences at the whole-population level, and the connection between these levels is mediated through adequate epidemiological models

Epidemiological models on the web

http://user.mendelu.cz/marik/dmb/sir.html

Some applications

Dynamics of HIV/AIDS

This example comes from Otto and Day (2007) and originally from Blower et al. (2000) AIDS results from the deterioration of the immune system and is the final stage of infection by HIV HIV is transmitted through the exchange of bodily fluids, predominantly through unprotected sexual intercourse, but also through sharing of unsterilized needles or transfusion with infected blood supplies UNAIDS & WHO: 33.2 million people worldwide were estimated to be infected by HIV in 2007, and 2.1 million people, including 330,000 children, were estimated to be killed by AIDS

Dynamics of HIV/AIDS Blower et al. (2000): model of spread of HIV among gay male community in San Francisco Motivation: concern that effective antiretroviral therapies (ART) might cause people to be less cautious when engaging in behavior posing a risk for HIV transmission Question: what this might mean for a further progression of the disease relative to when no ART was available Otto and Day (2007): simplified model version (ignores evolution of HIV resistance) Major assumption: mean number of sexual partners with whom an HIV– individual has unprotected sex per year increased from c before ART was available to c(1 + w ), for a constant w > 0

Dynamics of HIV/AIDS Three classes of gay males: uninfected individuals X (t), HIV+ individuals taking ART YT (t), and untreated HIV+ individuals YU (t) dX /dt = π − c(1 + w )λX − µX

dYU /dt = c(1 + w )λX + gYT − σYU − µYU − νU YU

dYT /dt = −gYT + σYU − µYT − νT YT

π . . . rate at which HIV– men enter the gay community in SF µ . . . rate at which gay men leave the community σ . . . rate at which untreated HIV+ men enter treatment g . . . rate at which treated HIV+ men abandon treatment νU (νT ) . . . death rate of un(treated) HIV+ men from AIDS

Dynamics of HIV/AIDS Force of infection λ (per capita rate at which susceptibles catch the disease) = probability that a sexual partner is HIV+ times the probability of acquiring HIV from this partner during sex, summed over the treated and untreated classes: λ = βU

YT YU + βT N N

N = X + YU + YT is the total number of gay males in the community Untreated partners are assumed to have a higher probability of HIV transmission than treated ones, βU > βT

Dynamics of HIV/AIDS 4

4

x 10

Gay community size

3.5 3 2.5 2 1.5 1 0.5 0 0

20

40 Time (years)

60

80

Dynamics of uninfected (red), treated HIV+ (green), and untreated HIV+ (blue) individuals when the risky behavior does not change (solid) and when it doubles (dash-dot); parameters as in Otto and Day (2007)

Dynamics of HIV/AIDS

Otto and Day (2007): 43% of the AIDS deaths that would have occurred within the gay population of San Francisco within the next ten years will be avoided by the use of ART if risky behavior does not increase (w = 0), but only 24% of deaths will be averted if risky behavior doubles (w = 1) From another perspective, there will be much higher probability of contracting HIV upon sexual intercourse when risky behavior doubles compared with when it does not increase

Hantavirus infection in rodents

This example comes from Allen et al. (2006) Hantavirus is primarily an infection in rodents that might be occasionally transmitted to humans Human infection occurs primarily through the inhalation of aerosolized saliva and/or excreta of infected rodents, or after individuals have been bitten by infected rodents Hemorrhagic fever with renal syndrome – HFRS (in Europe and Asia) or hantavirus pulmonary syndrome – HPS (in North and South America); mortality rate for HPS in the United States is 37%

Hantavirus infection in rodents

Need for two-sex models – two aspects of the disease call for an explicit consideration of the host sex structure: As males are more aggressive than females, contacts between two males generally result in greater transmission of the disease than contacts between two females or a male and a female (this implies that males have a higher disease prevalence than females) Infectious period is longer for males than for females Life expectancy of rodents is relatively short ⇒ latent period of the disease needs to be accounted for

Hantavirus infection in rodents Model classes: susceptible (S), exposed (E ), infectious (I ), and recovered (R), and then males (subscript m) and females (subscript f ) SEIR model for male rodents: dSm B(Nm , Nf ) = − βf Sm If + βm Sm Im − d(N)Sm dt 2 dEm = βf Sm If + βm Sm Im − σEm − d(N)Em dt dIm = σEm − γm Im − d(N)Im dt dRm = γm Im − d(N)Rm dt

Hantavirus infection in rodents SEIR model for female rodents: B(Nm , Nf ) dSf = − βf Sf If + βmf Sf Im − d(N)Sf dt 2 dEf = βf Sf If + βmf Sf Im − σEf − d(N)Ef dt dIf = σEf − γf If − d(N)If dt dRf = γf If − d(N)Rf dt Nm = Sm + Em + Im + Rm . . . total density of males Nf = Sf + Ef + If + Rf . . . total density of females N = Nm + Nf . . . total population density

Hantavirus infection in rodents

B(Nm , Nf ) . . . birth function Assuming 1:1 sex ratio of newborns, one of the most commonly used birth functions in sex-structured models is B(Nm , Nf ) = b

Nm Nf (Nm + Nf )/2

in which b is the mean litter size and Nm Nf /[(Nm + Nf )/2] is the pair formation rate Density-dependent death rate d(N) = d + d1 N, with 0 < d < b/2 and d1 > 0 is the same for males and females

Hantavirus infection in rodents

Requirements on parameters: Parameter βf scales the disease transmission rate between an infective female and a susceptible female or a susceptible male, βmf scales the disease transmission rate between an infective male and a susceptible female, and βm scales the disease transmission rate between an infective male and a susceptible male The above stated differences between male and female epidemiology translate to the following conditions: βm ≥ βmf ≥ βf and γf > γm

Hantavirus infection in rodents Disease-free equilibrium (DFE): Sm = K /2 = Sf and Em = Im = Rm = Ef = If = Rf = 0 No disease-induced mortality ⇒ every equilibrium must satisfy Sm + Em + Im + Rm = Sf + Ef + If + Rf = K /2 Basic reproduction number (using the approach by van den Driessche and Watmough 2002):  βm βf σK /4 + + R0 = b/2 + σ b/2 + γm b/2 + γf # p [βm (b/2 + γf ) + βf (b/2 + γm )]2 − 4βf (βm − βmf )(b/2 + γf )(b/2 + γm ) (b/2 + γm )(b/2 + γf )

Hantavirus infection in rodents

Importantly, R0 is proportional to the environmental carrying capacity of rodents K . As K increases, R0 also increases and so does the chance of a disease invasion. This relationship between R0 and K is a consequence of the assumption of mass action incidence, a reasonable assumption for rodent populations. This result agrees well with what happened in New Mexico in 1993. There, the outbreak of Sin Nombre virus (a hantavirus strain) was associated with increased densities of deer mice Peromyscus maniculatus, its primary host.

More complex microparasitic models

More complex microparasitic models

More realistic = more useful in applications Built upon simple models considered so far Model types we will consider: Models with free-living infectious stages Vector-host models Multi-host models Multi-strain models Multi-group models Stage-progression models

Models with free-living infectious stages

Models with free-living infectious stages

Some microparasites form free-living infectious stages Spores of bacteria, protozoa and fungi Capsules, polyhedra or free particles of viruses Important part of the parasite life cycle, carrying the infection from one host to the next Transmission now occurs between susceptibles and free-living infectious particles of the parasite Minimal model = equations for susceptibles and infectives, and for the free-living infectious stage of the parasite

Models with free-living infectious stages

Conceptual equation for the free-living infectious particles: Rate of change of the number or density of particles = (particle discharge rate) × (number of infectious hosts) – particle decay rate – particle consumption rate by the host Infectives release free-living infectious particles into the environment at a steady rate per infectious host Particles remain in the environment, with their numbers decaying away at a certain fixed rate Particles decrease in number also due to consumption by the hosts (both susceptible and infectious)

Models with free-living infectious stages

SI model with mass action incidence and exponential population growth (b > d, N = S + I ): dS = bN − βWS − dS dt dI = βWS − dI − αI dt dW = λI − µW − βWN dt W . . . size or density of free-living infectious particles in the environment

Models with free-living infectious stages

Models with free-living infectious stages

Reality: For many relevant parasites, Λ infectious particles are released only upon the death of an infective Modeling: A host releasing Λ infectious particles upon its death is essentially equivalent to a host releasing particles at a steady rate λ = Λ(α + d) throughout the expected lifespan 1/(α + d) of its infection

Models with free-living infectious stages

Anderson and May (1981)

Models with free-living infectious stages

Anderson and May (1981)

Models with free-living infectious stages R0 = (particle discharge rate) x (lifespan of infectious hosts) x (transmission rate per infectious particle) x (lifespan of particles) Number of free-living infectious particles released by an infectious host into the environment Number of susceptible hosts infected by any of these particles in the environment λ βN 1 1 R0 = λ βS = d +α µ + βN S=N d + α µ + βN

β influenced by numerous factors including host’s mobility and the chance that infection occurs after ingestion of (or other contact with) a free-living infectious particle

Models with free-living infectious stages If λ < d + α then R0 < 1 and the disease always dies out irrespectively of the host population size If λ > d + α then R0 = 1 implies threshold host density NT =

µ d +α β λ − (d + α)

and there are two possibilities:

if N > NT the disease is able to spread if N < NT the host population is too rare to support the disease, which is dying out; the host population will therefore grow exponentially at rate r = b − d until it exceeds NT , whereupon the disease is able to spread

Models with free-living infectious stages Four basic regimes of dynamical behavior: Parasites cannot persist and the host population grows exponentially at the disease-free rate r = b − d

Parasites are unable to regulate the host population but persist and are slowing down the rate of exponential host growth Disease regulates the host population to a stable equilibrium Disease regulates the host population to a stable limit cycle

This last possibility is a dynamical regime that we have not seen so far – regular host population fluctuations due to internal system dynamics

Models with free-living infectious stages

50

0 0

slower growth

λ

100

limit cycles stable

5 α

disease extinction 10

Highly pathogenic microparasites producing very large numbers of long-lived infectious stages are likely to cause non-seasonal cyclic changes in the abundance of their hosts and in the prevalence of infection

Models with free-living infectious stages

Population size or density

100

T

80 60 40 20 0 0

Disease prevalence I/N

120 total host free−living parasite critical total host N

5

10 Time

15

T

80 60 40 20

0.6

0.6

0.5

0.5

0.4 0.3 0.2 0.1 0 0

total host free−living parasite critical total host N

100

0 180

20

Disease prevalence I/N

Population size or density

120

185

190 Time

195

200

185

190 Time

195

200

0.4 0.3 0.2 0.1

5

10 Time

15

20

0 180

Models with free-living infectious stages

The peak in infection prevalence occurs shortly after the peak in host abundance The host population falls below the threshold value NT for a part of the cycle; the infection here survives primarily by virtue of its relatively long-lived transmission stages When N < NT , the infection prevalence declines effectively to zero, so that the disease seems to have disappeared form its host population even though the population of free-living infectious stages remains relatively abundant

Application

Larch budmoth Zeiraphera diniana and granulosis virus

Application

Anderson and May (1980): Without pathogens, moth population doubles in 8-9 months Without pathogens, moth life expectancy is about 3.5 months Pathogen average lifetime is 3 to 4 months Infected moths die on average 3 weeks after infection Pathogen discharge rate is 106 per year Transmission coefficient is set (arbitrarily) at 10−10 per year

Application

1

10 8 6 4

total host free−living parasite critical total host NT

2 0 600

700

800 Time (months)

900

1000

Disease prevalence I/N

Population size (log scale)

12

0.8 0.6 0.4 0.2 0 600

700

800 Time (months)

900

1000

The cyclic pattern of host abundance tends to be characterized by a slow rise and a rapid fall, whereas the cycle in the disease prevalence tends to be more symmetrical

Variation

Boots (1999): Rates of uptake of the infectious stages by susceptibles and infectives differ – reduced uptake rate for infectives dS = bN − βWS − dS dt dI = βWS − dI − αI dt dW = λI − µW − β(S + (1 − q)I )W dt

Variation

Four levels of debilitation: (a) q = 0, (b) q = 0.25, (c) q = 0.75, (d) q = 1

Vector-host models

Vector-host models

Fruitful area of research, as many important diseases of humans are transmitted by vectors and especially mosquitoes (e.g. malaria, bubonic plague, Lyme disease) New element: the need to include population dynamics of the vector Hosts infected by contacts with infectious vectors, and vectors in turn infected by contacts with infectious hosts Equations needed for each class of hosts and each class of vectors

Vector-host models

Simple model One host and one vector SIS model for hosts and SI model for vectors Four compartments: susceptible hosts (SH ), infectious hosts (IH ), susceptible vectors (SV ) and infectious vectors (IV ) Transmission rates assumed of the standard incidence type (vectors assumed to bite a fixed number of hosts per unit time)

Vector-host models

SH dSH = µNH − βV IV − µSH + γIH dt NH SH dIH = βV IV − µIH − γIH dt NH IH dSV = cNV − βH SV − cSV dt NH IH dIV = βH SV − cIV dt NH NH = SH + IH and NV = SV + IV are total population densities of hosts and vectors, respectively Both assumed constant here (births and natural deaths balanced and no disease-induced mortality)

Vector-host models Basic reproduction number (van den Driessche and Watmough 2002): s βH βV NV R0 = µ + γ c NH Near DFE (NH , 0, NV , 0): each infectious vector produces βV (SH /NH )/c|SH =NH = βV /c new infectious hosts over its mean infectious period 1/c each infectious host produces βH (SV /NH )/(µ + γ)|SV =NV new infectious vectors over its mean infectious period is 1/(µ + γ) the product gives the total basic reproduction number from vector to vector or from host to host

Vector-host models R0 < 1 if and only if βH βV NH > NV µ+γ c If the ratio of the total host population NH to the total vector population NV becomes large, the infection dies out This is because when there are many hosts (such as humans) relative to vectors (such as mosquitoes), the chance of someone being bitten twice in a quick succession – once to catch the infection and once to pass it on before recovery – is very small For a vector-borne infection to successfully invade the host population, the ratio of vectors to hosts has to be sufficiently large so that double bites are common

Vector-host models Analysis Rescaling sH = SH /NH , iH = IH /NH , sV = SV /NV , iV = IV /NV and setting βˆV = βV (NV /NH ): dsH /dt = µ − βˆV iV sH − µsH + γiH diH /dt = βˆV iV sH − µiH − γiH

dsV /dt = c − βH sV iH − csV diV /dt = βH sV iH − ciV Incorporating sH = 1 − iH and sV = 1 − iV :

diH /dt = βˆV iV (1 − iH ) − µiH − γiH diV /dt = βH (1 − iV )iH − ciV

Vector-host models Two equilibrium points: DFE (0, 0) → (NH , 0, NV , 0) Endemic equilibrium (iH∗ , iV∗ )

=

βH βˆV − c(b + γ) βH βˆV − c(b + γ) , βH (βˆV + b + γ) βˆV (βH + c)

!

SH∗ = (1 − iH∗ )NH , IH∗ = iH∗ NH , SV∗ = (1 − iV∗ )NV , IV∗ = iV∗ NV , βˆV = βV (NV /NH ) Local stability analysis R0 < 1: there is only the DFE (NH , 0, NV , 0) which is locally stable R0 > 1: the DFE is unstable and there is a unique endemic equilibrium (SH∗ , IH∗ , SV∗ , IV∗ ) which is locally stable

West Nile virus (WNV) infection This example comes from Wonham et al. (2004) Primary hosts of WNV are birds, especially crows, and the virus is transmitted via mosquitoes of the genus Culex Mammals (e.g. horses and humans) are secondary hosts, generally considered unimportant to disease persistence WNV is widespread in Africa, the Middle East and western Asia, with occasional European outbreaks introduced by migrating birds North America: the first recorded epidemic was initially detected in the New York state in 1999 and spread rapidly across the continent, with many bird, horse and human deaths left behind Control strategies focus primarily on the eradication of vector mosquitoes

West Nile virus (WNV) infection The simplest possible biologically relevant model for WNV transmission contains three compartments for birds and four compartments for mosquitoes

West Nile virus (WNV) infection

Birds dSB /dt = −ab IM SB /NB

dIB /dt = ab IM SB /NB − µV IB − gIB

dRB /dt = gIB Mosquitoes

dLM /dt = βM (SM + EM + IM ) − mLM − µL LM

dSM /dt = −ac SM IB /NB + mLM − µA SM

dEM /dt = ac SM IB /NB − kEM − µA EM

dIM /dt = kEM − µA IM

West Nile virus (WNV) infection Parameterization The American crow, Corvus brachyrhynchos, the bird that has suffered some of the highest mortality in the North American WNV epidemic, and the mosquito Culex pipiens sspp., a major North American WNV vector

West Nile virus (WNV) infection Basic reproduction number (using mean parameter values): s s k ab SM (0) ac SM (0) ≈ 0.465 R0 = µA NB (0) µV + g k + µA NB (0) When R0 < 1 the DFE is locally stable; when R0 > 1 it is locally unstable and an epidemic outbreak occurs The virus will invade if SM (0)/NB (0) > 4.625 For New York in 2000, a 40–70% reduction of the initial mosquito population, i.e. reducing SM (0)/NB (0) from 7.5–15 to less than about 4.6, would have prevented the WNV outbreak (bird control would have had the opposite effect)

Multi-host models

Multi-host models Many diseases are host-specific, but many others can infect multiple hosts Classic examples: foot-and-mouth disease, which can infect a large variety of livestock species such as cattle, sheep and pigs, and bovine tuberculosis, affecting badgers and cattle Different species may have different epidemiological responses to the same infection: in one species an infection may be short-lived and highly virulent, whereas in another species relatively long-term or mild infection can be the norm Multi-host infections pose a complex problem that can only be understood by means of mathematical models

Multi-host models Consider n host species and a single disease that can be transmitted both within and between the species SIS model (j = 1, . . . , n) n

X dSj Ik Sj = bj (Nj )Nj − − dj (Nj )Sj + γj Ij βjk (Nj ) dt Nj k=1

dIj = dt

n X k=1

βjk (Nj )

Ik Sj − dj (Nj )Ij − αj Ij − γj Ij Nj

In this model, species are assumed to interact only through the infection by a common pathogen – there are no competitive or predator-prey relationships

Multi-host models Many multi-host diseases have one host as a reservoir species and the other hosts as spillover species The highest transmission rate within the reservoir species, higher transmission rate from the reservoir species 1 than from the spillover species 2, . . . , n or between the spillover species:   βjk (Nj ), j, k = 2, . . . , n β11 (N1 ) > βj1 (Nj ) ≥  β (N ), j = 2, . . . , n 1j 1

The reservoir species has a longer period of infectivity than the spillover species: γ1 + α1 ≤ γj + αj , j = 2, . . . , n

Multi-host models Basic reproduction number (n)

R0

= spectral radius of the n × n matrix (Rjk )nj,k=1 where Rjk =

βjk (Kj ) , Kj = carrying capacity of species j γk + αk + bk

Rjj = basic reproduction number for species j, when there is no contact with the other species (j = 1, 2, . . . , n) It follows that R11 > Rjj for all j = 2, . . . , n Special case n = 2 (two host species):   q 1 R0 = (R11 + R22 ) + (R11 + R22 )2 + 4R12 R21 2

Multi-host models: Spillover species If a new spillover species is added to the system and all parameters of the original system with n hosts remain unchanged, then the (n+1) new system with n + 1 hosts has basic reproduction number R0 that satisfies (n+1)

R0

(n)

≥ R0

(1)

≥ R0

= R11

For a system with R0 < 1, addition of a new spillover species into the system may increase the reproduction number to a value greater than one, R0 > 1 Possibility of disease emergence increases with the number of host species ⇒ presence of spillover species can be an important factor in the emergence and persistence of infections

Multi-host models: Spillover species One reservoir host and one spillover species

McCormack and Allen (2007)

Multi-host models: Spillover species Another spillover species added

McCormack and Allen (2007)

Multi-host models: Dilution effect Tick-borne virus and two host species (Norman et al. 1999) One host non-viraemic (suffers no effects from the disease) and the other host viraemic (suffers mortality due to the disease) Ticks need to feed on hosts to grow to the next developmental stage and to reproduce (adult females) Tick developmental stages: larvae (L), nymphs (N), adults (A) Ticks born as susceptibles; nymphs and adults can be either susceptible (Ns , As ) or infectious (Ni , Ai ) Non-viraemic hosts kept constant (H1 ), viraemic hosts classified as either susceptible (H2s ), infectious (H2i ), or recovered (H2r ) Both ticks and viraemic hosts subject to logistic growth

Multi-host models: Dilution effect

dL/dt = (β5 D5 H1 + β6 D6 H2 )(Ai + As )(a1 − s1 T ) − bL − β1 H1 L − β2 H2 L

dNi /dt = β2 D2 H2i L − bNi − β4 H2 Ni − β3 H1 Ni

dNs /dt = β2 D2 (H2s + H2r )L + β1 D1 H1 L − bNs − β4 H2 Ns − β3 H1 Ns

dAi /dt = β4 D4 H2 Ni + β4 D4 H2i Ns + β3 D3 H1 Ni − bAi − β5 H1 Ai − β6 H2 Ai dAs /dt = β4 D4 (H2s + H2r )Ns + β3 D3 H1 Ns − bAs − β5 H1 As − β6 H2 As

dH2s /dt = (a2 − s2 H2 )H2 − b2 H2s − β7 H2s Ni − β8 H2s Ai

dH2i /dt = β7 H2s Ni − β8 H2s Ai − (b2 + α + γ)H2i

dH2r /dt = γH2i − b2 H2r

β’s . . . probabilities of tick stages feeding on hosts D’s . . . proportions of ticks that feed and moult successfully T = L + Ns + Ni + As + Ai , H2 = H2s + H2i + H2r

Multi-host models: Dilution effect

The non-viraemic hosts can have significant impact on the viraemic host. Either they amplify the tick population and cause the virus to persist, or they dilute the infection and cause it to die out.

Invasion of gray squirrels across Europe

System with the red versus gray squirrel and parapoxvirus (Tompkins et al. 2003)

Invasion of gray squirrels across Europe

Gray squirrels were introduced from America to Europe at the beginning of the 20th century They have displaced red squirrels from much of its native range Gray squirrels have an innate competitive advantage over red ones But this advantage is not sufficient to explain the gray’s rapid expansion and the accompanying red’s decline Parapoxvirus was postulated to speed up competitive replacement A two-host epidemiological model is required to understand and predict a likely competitive outcome

Invasion of gray squirrels across Europe SIR (for gray) – SI (for red) model with interspecific competition   N G + cR N R dSG NG − bG SG − (βGG IG + βGR IR )SG = aG − dt KG dIG = (βGG IG + βGR IR )SG − bG IG − γG IG dt dRG = γG IG − bG RG dt   dSR N R + cG N G NR − bR SR − (βRG IG + βRR IR )SR = aR − dt KR dIR = (βRG IG + βRR IR )SR − bR IR − αR IR dt

Spatial scale: an individual forest

Invasion of gray squirrels across Europe

4 3.5 Squirrel density

3

red without parapox gray without parapox red with parapox gray with parapox

2.5 2 1.5 1 0.5 0 0

5

10 Time (years)

15

20

Parameters: Tompkins et al. (2003) Even the simple, non-spatial and deterministic model highlights the importance of parasites in affecting competition

Invasion of gray squirrels across Europe More realistic model: (1) spatial – invasion across a country, (2) stochastic – invasions and extinctions involve a low number of individuals

Tompkins et al. (2003) (a) observations (b) competition and infection (c) only competition

Multi-strain models

Multi-strain models Many infections that we consider as a single disease are in reality comprised of multiple strains, which interact through cross-immunity (or not) that is invoked within the host – malaria, dengue fever, influenza, etc. Multi-strain models = models of a host population that can be (simultaneously) infected by multiple strains of a disease Assumptions about conferred immunity and interaction between pathogen strains: either transmission of or susceptibility to one strain is modified due to passed infection by another strain Need to uniquely identify the entire infection history of individuals and their immunological status with respect to the various pathogen strains under consideration

Multi-strain models

Strains compete for limited resources = susceptible hosts Application Genetically modified (GM) viruses used to control some pest populations – GM viruses compete for susceptible hosts with naturally occurring viruses used to create the GM ones Virus-vectored immunocontraception (VVIC): Used to control mammal pest populations A naturally occurring virus is modified to sterilize one or the other sex of the pest species Example: myxomatosis to control rabbits in Australia

Multi-strain models

Types of TWO-strain models we will consider: Complete cross-immunity Partial cross-immunity No cross-immunity Superinfection Enhanced susceptibility Major question: can various diseases or disease strains coexist in the host?

Complete cross-immunity (Life-long) immunity to one strain automatically confers (life-long) immunity to the other Two classes of infectives, each corresponding to one strain SIR model dS/dt = bN − β1 SI1 /N − β2 SI2 /N − dS dI1 /dt = β1 SI1 /N − γ1 I1 − dI1

dI2 /dt = β2 SI2 /N − γ2 I2 − dI2

dR/dt = γ1 I1 + γ2 I2 − dR

b = d ⇒ constant total host population size N = S + I1 + I2 + R

Complete cross-immunity Strain-specific basic reproduction number (k)

R0

= βk /(d + γk ), k = 1, 2

Equations for I1 and I2 give conflicting conditions for equilibrium S ⇒ there is no strain coexistence equilibrium Numerical simulations ⇒ at most one strain can persist in this model, coexistence is not possible When competing strains provide complete protection for each (k) (k) other and R0 > 0 for both strains, the strain with the larger R0 will force the other strain to extinction

Complete cross-immunity (k)

Although R0 determines the eventual competitive outcome, for S(0) ∼ N, strain with a larger per capita growth rate βk − d − γk (k) will be initially favored, even if it has lower R0 and is eliminated

β1 /β2 > (d + γ1 )/(d + γ2 ) and

β1 − β2 < (d + γ1 ) − (d + γ2 )

Strains with different transmission modes Thrall and Antonovics (1997) Two strains, one spreading under standard incidence (and causing some sterilization) and the other under mass action incidence (and causing disease-induced mortality) SI model dS/dt = b[S + I1 + (1 − )I2 ] − β1 SI1 − β2 SI2 /N − dS

dI1 /dt = β1 SI1 − dI1 − αI1 dI2 /dt = β2 SI2 /N − dI2

α . . . mortality induced by parasite strain 1  . . . degree of sterility induced by parasite strain 2 b = b0 /(1 + ξN) . . . negative density dependence in reproduction

Complete cross-immunity

grey = STD wins black = OID wins white = coexistence

Either mode can win, and there is also an area in the parameter space in which both transmission modes coexist; this is most likely when the standard incidence strain is mildly sterilizing and the mass action incidence strain causes intermediate levels of mortality

Partial cross-immunity Having experienced and recovered from one strain provides some form of limited protection against other strains Most commonly studied and most widely applicable form of immunity Range of differing assumptions can be made about the nature of partial cross-immunity: Protection can operate either through reduced susceptibility, reduced transmissibility, or a mixture of the two Protection might be conferred only to a fraction of recovered individuals (in what follows, we assume that all recovered individuals experience the same reduction)

Partial cross-immunity

Eight compartments, and a variety of different parameters which measure the type and extent of cross-immunity Notation: NSI – number of individuals that are susceptible to strain 1 and infected with strain 2 (similarly for the other combinations)

Partial cross-immunity SIR model (standard incidence, b = d) dNSS /dt = bN − dNSS − β1 NSS (NIS + a1 NIR )/N − β2 NSS (NSI + a2 NRI )/N

dNIS /dt = −dNIS + β1 NSS (NIS + a1 NIR )/N − γ1 NIS

dNRS /dt = −dNRS + γ1 NIS − α2 β2 NRS (NSI + a2 NRI )/N dNSI /dt = −dNSI + β2 NSS (NSI + a2 NRI )/N − γ2 NSI

dNRI /dt = −dNRI + α2 β2 NRS (NSI + a2 NRI )/N − γ2 NRI

dNSR /dt = −dNSR + γ2 NSI − α1 β1 NSR (NIS + a1 NIR )/N dNIR /dt = −dNIR + α1 β1 NSR (NIS + a1 NIR )/N − γ1 NIR

dNRR /dt = −dNRR + γ1 NIR + γ2 NRI

For individuals who have recovered from the other strain: αk . . . proportional reduction in susceptibility to strain k ak . . . proportional reduction in transmissibility of strain k

Partial cross-immunity Assumption: strain 1 successfully invaded the host (1) (R0 = β1 /(γ1 + d) > 1) and the system host-strain 1 is at the endemic equilibrium: ∗ NSS

γ1 + d ∗ , NIS =N =N β1



d d − γ1 + d β1



,

∗ NRS

=N



γ1 γ1 − γ1 + d β1



Denoting I2 = NSI + a2 NRI , strain 2 invades if dI2 /dt > 0 at that endemic equilibrium: dI2 = β2 I2 dt

"

1 (1)

R0



1 (2)

R0

!

#  γ1  (1) R0 − 1 > 0 + a2 α2 β1

Partial cross-immunity

Strain 2 thus invades if: (2)

R0

(2)

if R0

(1)

> R0

(1)

< R0

then it invades if: a2 α2 >

1 (1)

R0 − 1

1 (2)

R0



1 (1)

R0

!

that is, if the effects of cross-immunity are relatively weak

Partial cross-immunity Symmetry: we get analogous conditions by first allowing strain 2 to successfully invade the host and reach the host-strain 2 endemic equilibrium and then seeking when strain 1 can invade Coexistence of strains is possible when both can invade when the other is at the endemic equilibrium with the host (1)

(2)

Assumptions: R0 > 1 and R0 > 1 (both need to be at the (1) (2) endemic equilibrium with the host), R0 > R0 (strain 1 invades the host-strain 2 system) Coexistence of strains is possible (strain 2 invades the host-strain 1 (1) (2) system) when e.g. R0 ≈ R0 and cross-immunity is relatively weak (a2 α2 ≈ 1)

Partial cross-immunity

The invasibility condition contains only the product ak αk , so it does not matter for a two-strain model whether partial cross-immunity acts on susceptibility (α), transmissibility (a), or a mixture of the two Generally not true when more than two strains are considered When many disease strains are circulating within the population, different assumptions about the nature of cross-immunity can have profound effects on the model outcome Detailed understanding of the host’s immunological responses is necessary to deal with realistic multi-strain models

No cross-immunity When coinfection is assumed not to occur (often individuals infected by one infection are not mixing with enough individuals to catch subsequent infections): SIR model (standard incidence, b = d) dNSS /dt = bN − dNSS − β1 NSS (NIS + NIR )/N − β2 NSS (NSI + NRI )/N

dNIS /dt = −dNIS + β1 NSS (NIS + NIR )/N − γ1 NIS

dNRS /dt = −dNRS + γ1 NIS − β2 NRS (NSI + NRI )/N dNSI /dt = −dNSI + β2 NSS (NSI + NRI )/N − γ2 NSI

dNRI /dt = −dNRI + β2 NRS (NSI + NRI )/N − γ2 NRI

dNSR /dt = −dNSR + γ2 NSI − β1 NSR (NIS + NIR )/N

dNIR /dt = −dNIR + β1 NSR (NIS + NIR )/N − γ1 NIR

dNRR /dt = −dNRR + γ1 NIR + γ2 NRI

No cross-immunity

Essentially the model for partial cross-immunity, with ak = αk = 1 Condition derived for the case of partial cross-immunity can be used to assess whether coxistence is possible for no cross-immunity For strain 2 invading the host-strain 1 endemic equilibrium, we (2) (1) must have either R0 > R0 or otherwise ! 1 1 1 − (1) < 1 (1) (2) R0 − 1 R0 R0

No cross-immunity

When coinfection is allowed to occur: we have two independent infections that do not interact with one another Coexistence of two strains is possible whenever both strain-specific (1) basic reproduction numbers allow for strain invasion: R0 > 1 and (2) R0 > 1

Superinfection One strain may “superinfect” an individual already infected with another strain SIS model to describe each strain in absence of the other

Superinfection gives rise to new infection in the corresponding infectious compartment

Superinfection

Let strain 1 superinfect an individual already infected with strain 2, giving rise to new infection in compartment I1 Let ν > 0 be the transmission rate for superinfection SIS model (standard incidence, b = d) dS/dt = bN − dS + γ1 I1 + γ2 I2 − β1 SI1 /N − β2 SI2 /N

dI1 /dt = β1 SI1 /N − dI1 − γ1 I1 + νI1 I2 /N

dI2 /dt = β2 SI2 /N − dI2 − γ2 I2 − νI1 I2 /N

ν = 0 . . . complete cross-immunity ⇒ no strain coexistence

Superinfection

Strain-specific basic reproduction number: (k)

R0

= βk /(d + γk ), k = 1, 2

The host-strain 2 endemic equilibrium (2)

(2)

S ∗ /N = 1/R0 , I2∗ /N = 1 − 1/R0 Strain 1 can invade this equilibrium if 1 dI1 1 = β1 (2) − (d + γ1 ) + ν I1 dt S ∗ , I ∗ R 2

0

1−

1 (2)

R0

!

>0

Superinfection

This holds if and only if (1)

R12 ≡

R0

(2)

Interesting case: R0

(2)

R0

ν + d + γ1

1−

1 (2)

R0

!

>1

(1)

> 1 > R0 , but R12 > 1

Strain 2 can invade the fully susceptible host population when alone but strain 1 cannot, and yet strain 1 can invade the host-strain 2 endemic equilibrium Coexistence can in this case occur if ν is sufficiently large

Enhanced susceptibility

Presence of one infection in a host individual can increase susceptibility of this individual to other infections Classic example: sexually transmitted infections Sexually transmitted infections usually conform to an SIS model, where after treatment infectious individuals are again susceptible

Enhanced susceptibility SIS model (standard incidence, no demography) NSS (NIS + NII ) NSS (NSI + NII ) dNSS = −β1 − β2 + γ1 NIS + γ2 NSI + γ3 NII dt N N dNIS NIS (NSI + NII ) NSS (NIS + NII ) = β1 − γ1 NIS − βˆ2 dt N N NSS (NSI + NII ) NSI (NIS + NII ) dNSI = β2 − γ2 NSI − βˆ1 dt N N dNII N (N + N ) N (N + N ) SI IS II IS SI II = βˆ1 + βˆ2 − γ3 NII dt N N

Assumption: infections passed on independently (individuals with both infections do not pass both on to an individual they infect) Enhanced susceptibility: βˆ1 > β1 and βˆ2 > β2

Enhanced susceptibility (k)

Strain-specific basic reproduction number: R0

= βk /γk , k = 1, 2

At invasion, prevalence is low and coinfection is very rare Terms with βˆ1 , βˆ2 and γ3 do not play a role When R0 < 1 for both infections neither can invade a totally susceptible population Once both infections are established and coinfection is common Terms with βˆ1 , βˆ2 and γ3 can play a pivotal role in maintenance of both infections If βˆ1 and βˆ2 are relatively large, given a substantial prevalence of infection 1, the average susceptibility to infection 2 is increased and infection 2 can persist, and vice versa Infections cannot invade but may persist once established

Enhanced susceptibility

Enhanced susceptibility Implications for disease control It may be erroneous to study sexually transmitted infections in isolation and a multi-strain approach may be necessary Small changes in βˆ1 and βˆ2 , i.e. in sexual practices, may cause significant changes in the prevalence Presence of a sexually transmitted infection may make it difficult to eradicate the others compared with a combined control strategy On the other hand, reducing the prevalence of one infection may lead to a reduction of the prevalence of the other

Multi-strain models and evolution

Multi-strain models are of high importance from an evolutionary perspective, for studies of disease evolution Multi-strain models = two or more strains compete for survival in a common host Evolution = a rare mutant strain (more or less virulent, more resistant to treatment, etc.) attempts to invade a resident strain established at an endemic equilibrium with its host The presented models do not comprise explicit genetics; genetics is only implicit in, for example, the differing rates of transmission or virulence of different virus strains

Multi-group models

Multi-group models Division of heterogeneous host population into several homogeneous groups based on individual behavior Each group subdivided into common epidemiological classes Majority of multi-group models are built for sexually transmitted diseases, such as HIV/AIDS or gonorrhea, where behavior is an important factor in the probability of contracting the disease, and where the behaviorally distinguished groups represent the sexual risk groups such as prostitutes, homosexuals or celibate individuals Models of this kind are sometimes also referred to as risk-structured models, given that different groups have different risk of contracting the disease

Multi-group models

Two-group SIS model with high-risk and low-risk individuals (standard incidence, no demography) SH IH SH IL dSH = −βHH − βHL + γH IH dt N N dIH SH IH SH IL = βHH + βHL − γH IH dt N N SL IH SL IL dSL = −βLH − βLL + γL IL dt N N dIL SL IH SL IL = βLH + βLL − γL IL dt N N

Constant number of individuals in each group: N = NH + NL

Multi-group models

Four transmission parameters, representing mixing between different (risk, social or sexual) groups These can be grouped into the WAIFW (Who Acquires Infection From Whom) matrix:   βHH βHL β= βLH βLL This matrix plays a role analogous to that the scalar parameter β plays in unstructured (with respect to groups) epidemiological models

Multi-group models Assumptions on β Individuals in the high-risk group are at a higher risk of being infected: βHH + βHL > βLH + βLL Assortative mixing: individuals from the high-risk group are more likely to have sex with other high-risk individuals, and analogously for low-risk individuals – diagonal elements of the WAIFW matrix dominate, βHH , βLL > βHL , βLH , with βHH > βLL Interactions between groups are symmetric: βHL = βLH (this does not need to hold generally, for example, when one group is, for some reason, more susceptible)

Multi-group models Example 

 NL 10 0.1 NH = 0.2, = 0.8, γ = 1 β= , 0.1 1 N N 0.12 Log proportion infectives

Proportion infectives

0.1

0 high−risk low−risk

0.08 0.06 0.04 0.02 0 0

5

10 Time

15

20

−1

high−risk low−risk

−2 −3 −4 −5 −6 0

5

10 Time

15

Multi-group models

Three phases 1

Transient phase: dependence on initial conditions

2

Exponential growth phase: driven by R0

3

Equilibrium phase: approach to a stable endemic equilibrium

Multi-group models Basic reproduction number R0 for individual reproductive groups: plays a role when infection starts in a single group (H)

High-risk group: R0

(L)

Low-risk group: R0

= (βHH NH + βHL NH )/(NγH ) = 2.02

= (βLH NL + βLL NL )/(NγL ) = 0.88

R0 for the entire system is the dominant eigenvalue of the matrix   βHH NH /(NγH ) βHL NH /(NγL ) G= which is 2.0013 βLH NL /(NγH ) βLL NL /(NγL ) R0 for the entire population is bounded by values for each group Exponential growth rate = R0 − 1

Multi-group models: Application to HIV Simplified dynamics of HIV among male homosexuals (Keeling and Rohani 2008, pages 71-74) Five groups with 0, 3, 10, 60, and 100 partners per year in each group, with sizes (0.06 0.31 0.52 0.08 0.03)N and transmission parameters  0  0  β = 0.0217   0  0 0

0 0.65 2.15 12.9 21.5

 0 0 0 2.15 12.9 21.5   7.17 43.1 71.8   43.1 258 431  71.8 431 718

Multi-group models: Application to HIV Equations (k = 1, . . . , 5, b = d) 5

X Sk Ij dSk βkj = bNk − − dSk dt N j=1

5 X

Sk Ij dIk βkj = − dIk − γIk dt N j=1 dAk = γµIk − dAk − αAk dt

A . . . group of infected individuals that have developed AIDS µ . . . proportion of infected individuals who develop AIDS Individuals with AIDS are assumed to inhibit their sexual behavior and therefore no longer contribute to the infection spread

Multi-group models: Application to HIV

Other parameters: d = 0.0312, γ = 0.2, α = 1, µ = 0.3

Stage progression models

Stage progression models

Class I divided into a sequence of subclasses through which the disease progresses (can be applied to any other relevant class such as E ) This subdivision allows to account for diseases such as HIV/AIDS in which the risk of transmission varies with time since infection

Stage progression models SI model (standard incidence, b = d = dk ) m−1 X SIk dS βk = bN − dS − dt N k=1

m−1 X SIk dI1 βk = − ν1 I1 − d1 I1 dt N k=1

···

dIk = νk−1 Ik−1 − νk Ik − dk Ik , k = 2, . . . , m − 1 dt dIm = νm−1 Im−1 − dm Im dt

The final stage m assumed to have zero infectivity due to morbidity

Stage progression models

R0 =

β1 β2 ν1 βm−1 ν1 . . . νm−2 + + ··· + ν1 + d1 (ν1 + d1 )(ν2 + d2 ) (ν1 + d1 ) . . . (νm−1 + dm−1 )

kth term – number of new infections produced by a typical individual during the time it spends in the kth infectious stage νk−1 /(νk−1 + dk−1 ) . . . fraction of individuals reaching stage k − 1 that progress to stage k 1/(νk + dk ) . . . average time an individual entering stage k spends in stage k kth term – product of the infectivity of individuals in stage k, fraction of initially infected individuals surviving at least to stage k, and average infectious period of an individual in stage k

Modeling recovery

Standard assumption: recovery rate scales linearly with the number or density of infectives: rate of recovery = γI γ is the (per capita) recovery rate What does this mean from the perspective of an infectious individual? What is the mean length of the infectious period? Consider a cohort of hosts who are all infectious at time 0, prevent further infections, and let I (τ ) denote the number of hosts who are still infectious at time τ > 0

Modeling recovery With the linear recovery rate, the number of infectives changes as dI /dτ = −γI Solution I (τ ) = I (0)e −γτ ⇒ e −γτ is the proportion of infectives remaining infective at time τ The probability that an individual recovers before time τ after contracting the disease thus is P(recovery time < τ ) = 1 − e −γτ The mean length of the infectious period of an individual: Z ∞ 1 τ γe −γτ dτ = γ 0

Modeling recovery

The linear recovery rate γI thus means that the mean time an individual is infectious is 1/γ and the individual times of infection are exponentially distributed around this mean Not much realistic but mathematically “elegant” Although the standard class-structured models, with their associated exponential distributions, are a useful tool for understanding the dynamics of infection, when very precise predictions are required it becomes necessary to include much more of the available information; details often matter How to include more realistic distributions of the infectious period?

Modeling recovery

Real data: it seems that the infectious period does not follow an exponential distribution but rather seems to be of constant duration To account for such more realistic distributions, we need to relax the assumption that the probability of recovery does not depend on the time since infection There are several ways to do that, of which the simplest one is the method of stages The basic idea of the method of stages is to replace the infective compartment by a series of n successive infective compartments, each with an exponential distribution of the same parameter γn

Modeling recovery The total duration of the infectious period = the sum of n identical and independent exponential distributions = a gamma distribution of the infectious durations with probability density function f (τ ) =

(γn)n n−1 −γnτ τ e Γ(n)

where Γ(n) is the gamma function n = 1 . . . the exponential distribution n large . . . the gamma distribution tends toward a normal one n → ∞ . . . the delta (fixed duration with no variance) distribution Scaling from exponential (when there is just one subclass) to relatively constant (when there are many subclasses)

Modeling recovery

Gamma distribution for γ = 0.5 and various values of the number n of classes. For all values of n the mean is equal to 1/γ = 2.

Modeling recovery: Application to SARS Keeling and Rohani (2008), pages 100-102

Four model classes: susceptible, exposed (3 subclasses), infectious (3 subclasses), hospitalized (10 subclasses)

Seasonality

Seasonality Childhood infections

Greater mixing of children during school terms than school holidays

Seasonality

Many diseases are affected by seasonal variations Outbreaks of influenza most commonly occur in winter Disease transmission may vary with temperature . . . or with other climatic factors such as moisture which affects e.g. biting insects . . . or with changes in host behavior throughout the year: increased contacts arising from flocking behavior, seasonal migration, and breeding aggregations Diseases of animals may also be influenced by patterns of seasonal birth

Seasonality

In all these cases, we are dealing with a kind of forcing external to the infection This forcing can be conveniently modeled by letting a model parameter fluctuate in time General question: How do seasonally varying model parameters act as a forcing mechanism to affect infection dynamics? Specific question: What are the effects of varying amplitude of oscillations in model parameters on host population dynamics?

Climatic variations Easier spread of a disease in a portion of the year than in another one: transmission coefficient β is a function of time t: β(t) = β0 (1 + β1 cos(2πt/T ))

β0 . . . baseline transmission coefficient T . . . period of forcing (β(0) = β(T )), usually T = 1 year β1 . . . amplitude of seasonality relative to β0

Climatic variations

Endemic SIR model SI dS = bN − β(t) − dS dt N SI dI = β(t) − γI − dI dt N dR = γI − dR dt Assumption: b = d ⇒ N = S + I + R does not vary with time

Climatic variations Endemic SIR model Absence of seasonal forcing (β(t) = β0 ): there is a unique endemic equilibrium which is stable and frequently attained through damped oscillations

Size or density of infectives

2

1.5

1

0.5

0 0

10

20

30 40 Time (years)

50

60

Climatic variations Effects of varying the amplitude β1 on model behavior β1 = 0.08 biennial cycles

β1 = 0.05 annual cycles −4

β1 = 0.1 chaotic behavior −3

−4

x 10

x 10

5

14

x 10

3.4 12

3 2.8 2.6 2.4

10 8 6

2.2

4

2

2

1.8 190

192

194 196 Time (years)

198

200

4 Fraction infectives

Fraction infectives

Fraction infectives

3.2

190

3 2 1

192

194 196 Time (years)

198

200

190

192

194 196 Time (years)

198

200

Now we plot a bifurcation diagram: we record a system state, here infection prevalence, at a fixed (always the same) time instant within each cycle of period T

Climatic variations Effects of varying the amplitude β1 on model behavior

Annual proportion infectives

1

10−5

10−10

10−15 0

0.02

0.04 0.06 0.08 0.1 Amplitude of seasonality β1

0.12

one dot = annual cycles, two dots = biennial cycles, many dots = chaos

Dynamical complexity increases as seasonality becomes stronger: from annual epidemics to biennial outbreaks to chaos

Climatic variations Single bifurcation diagram does not need to fully capture the entire range of model behavior Multiple stable attractors, if any, imply that different initial conditions can lead to different bifurcation diagrams extrapolated initial conditions

fixed initial conditions 1 Annual proportion infectives

Annual proportion infectives

1

10−5

10−10

10−15 0

0.02

0.04 0.06 0.08 0.1 Amplitude of seasonality β1

0.12

10−5

10−10

10−15 0

0.02

0.04 0.06 0.08 0.1 Amplitude of seasonality β1

0.12

Climatic variations Why do some amplitudes of forcing result in annual cycles and some in multiennial ones? For the disease to invade the host population, we require   dI β(t) S SI = β(t) − γI − dI = I (γ + d) −1 >0 dt N γ+d N This is satisfied if s ≡ S/N > sT (t) ≡ (γ + d)/β(t) β varies with time ⇒ threshold fraction of susceptibles sT also varies with time This is the core of oscillatory dynamics: we may have s < sT (t) in a portion of the year and s > sT (t) in another one

Climatic variations

transmission coefficient β(t)

annual cycles

biennial cycles

Legend: thin line = sT ; thick line = s (dark for s > sT and light for s < sT ); dashed line = i (= I /N)

Climatic variations

As the amplitude of seasonality β1 increases, larger epidemics are generated that lower the fraction of susceptibles s such that recovery to above the threshold sT takes far longer, resulting in longer-period cycles

Childhood diseases School terms versus shool holidays Seasonally forced models of childhood diseases, such as measles, now more often assume β(t) is piece-wise constant, attempting to capture aggregation of children in schools: β(t) = β0 (1 + β1 Term(t))

Transmission rate, β(t)

1350

Term(t) = 1 during school terms (high transmission)

1300 1250

Term(t) = −1 during shool holidays (low transmission)

1200 1150 0

0.5

1 Time (years)

1.5

2

Childhood diseases

“Switching” During term-times β remains constant at the higher value and the trajectories approach the stable endemic equilibrium given by this transmission rate in exactly the same manner as predicted by an unforced model When holidays start, a new, lower rate β is operative and a new stable endemic equilibrium exists, with the trajectories heading towards it We may view the changing values of β in term-time-forced models as switching the model between two attracting equilibria

Childhood diseases

Term-time forcing in general produces a much simpler bifurcation diagram than sinusoidal forcing: Biennial outbreaks occur for a larger region of parameter space Chaotic dynamics only occur for relatively high values of β1

School terms versus shool holidays

Is the precise shape of seasonality important? Choice of functional form used to represent seasonality in the transmission term can have a substantial qualitative, as well as quantitative, effect on infectious disease dynamics

Wildlife populations: transmission rate

Seasonal changes in flocking and social mixing = periods of high transmission rates: Autumn and winter aggregations at bird feeders Breeding period of seals when animals aggregate on beaches (observed to coincide with outbreaks of phocine distemper virus in these animals) Seasonal migration in insects, birds and mammals Such situations can be dealt with in a similar manner to the seasonal changes experienced for human diseases (aggregation of hosts)

Wildlife populations: birth rate

Seasonality arises from concentrating host births into a short period of time “Pulse” of hosts recruited into the population at approximately the same time each year generates a “pulse” of susceptibles that enter the host population A growing number of empirical studies point to seasonal births, possibly in combination with changes in social behavior, as a factor important to the dynamics of wildlife pathogens, with examples including cowpox virus infecting voles and mice and phocine distemper virus infecting seals

Wildlife populations: birth rate Equation for susceptibles now becomes dS = b(t)N − βSI − dS dt Births occur in an annual pulse ( b0 T < t < T + L b(t) = 0 T +L < t 1, R02 < 1 disease starts in patch 1

disease starts in patch 2

30

30 S

S

20

1

25

I1 S2 I2

15 10 5 0 0

1

I

1

Population size

Population size

25

20

S2 I2

15 10 5

10

20

30 Time

40

50

60

0 0

10

20

30 Time

40

50

60

In general rate of disease spread depends on the starting patch

“Animal” metapopulation models SI model (k = 1, . . . , n): X X Ik dSk = bk Nk − βk Sk − dk Sk + mkj Sj − mjk Sk dt Nk j6=k j6=k X X Ik dIk = βk Sk − dk Ik − αk Ik + mkj Ij − mjk Ik dt Nk j6=k

j6=k

Coupling driven by parameters mkj = rate at which host migrate from patch j to k It is the direct movement of infectives that spreads the pathogen, since transmission between subpopulations occurs only following the movement of an animal

“Animal” metapopulation models Virus transmission on transports (Cui et al. 2006) How does transport-related infection affect dynamics of an infectious disease between two patches Standard studies – large country with a small number of potentially large cities and very sparse or even nonexistent rural population and good transportation system – movement from one city to another is fast, and the propagation of an epidemic takes place only at the destination location Some developing countries – dense crowd on trains, airplanes, and also sanitary conditions relatively bad – virus transmission on transports may be an important factor for the spread and persistence of infectious diseases

“Animal” metapopulation models

Assumptions m . . . rate at which hosts travel between patches ν . . . transmission efficiency during travel . . . when individuals travel from city j to city k, disease is assumed transmitted via standard incidence: ν

Sj Ij Sj Ij (αSj )(αIj ) = να = να αSj + αIj Sj + Ij Nj

= standard incidence with transmission coefficient να Individuals who are travelling do not give birth and do not die

“Animal” metapopulation models

SI model S1 I1 S2 I2 dS1 = bN1 − β − dS1 − αS1 + αS2 − να dt N1 N2 S1 I1 S2 I2 dI1 =β − dI1 − αI1 + αI2 + να dt N1 N2 dS2 S2 I2 S1 I1 = bN2 − β − dS2 − αS2 + αS1 − να dt N2 N1 dI2 S2 I2 S1 I1 =β − dI2 − αI2 + αI1 + να dt N2 N1

“Animal” metapopulation models

3000

1 * ν * I ν

2500

Disease prevalence, I*ν/(S*ν+I*ν)

Population size or density

S

2000 1500 1000 500 0 0

0.5 1 1.5 Transmission on transport, ν

2

0.8 0.6 0.4 0.2 0 0

0.5 1 1.5 Transmission on transport, ν

2

Infection on transports facilitates the infection spread and negatively affects the host population at endemic equilibrium

“Human” metapopulation models

Commuters Skj . . . number of susceptibles currently in city k that live in city j (analogously for members of the other epidemiological classes) ckj . . . rate at which individuals leave their home city j and commute to city k rkj . . . rate of return from city k of those that live in city j νkj . . . refers to individuals who live in city j but are born in city k P j Nkj . . . number of individuals currently in city k

“Human” metapopulation models SI model (k = 1, . . . , n): P X X dSkk j Ikj − dkk Skk − cjk Skk + rjk Sjk = νkk − βk Skk P dt j Nkj j j P dSkj j Ikj = νkj − βk Skj P − dkj Skj + ckj Sjj − rkj Sij , k 6= j dt j Nkj P X X dIkk j Ikj − dkk Ikk − cjk Ikk + rjk Ijk = βk Skk P dt j Nkj j j P dIkj j Ikj = βk Skj P − dkj Ikj + ckj Ijj − rkj Ikj , k 6= j dt j Nkj

In many situations the parameters c and r can be found from travel statistics

Influenza – Colizza et al. (2007b)

Application Impact of airline travel on the geographic spread of influenza

Influenza – Colizza et al. (2007b)

Influenza – Colizza et al. (2007b)

Sk Ikt Sk Iknt Sk Ika − β − r β β Nk Nk Nk + movement nt t Sk Ik S Ia S I + rβ β Nk kk − Lk + movement dLk /dt = β Nk kk + β N k dIkt /dt = pt (1 − pa )Lk − µIkt + movement dIknt /dt = (1 − pt )(1 − pa )Lk − µIknt dIka /dt = pa Lk − µIka + movement dRk /dt = µ(Ikt + Iknt + Ika ) + movement

dSk /dt = −β

β is assumed to fluctuate seasonally – easier transmission in winter, more difficult in summer – depends on the hemisphere

Influenza – Colizza et al. (2007b) IATA air travel data 3100 airports (urban areas) in 220 countries, 17182 flight connections = 99% worldwide traffic

Influenza – Colizza et al. (2007b)

Movement The number of passengers of category X traveling from a city j to a city k is an integer random variable, in that each of the Xj potential travellers has a probability pjk = wjk ∆t/Nj to go from j to k in the time interval ∆t wjk . . . traffic flow between cities j and k per unit time (IATA) The occupancy rate of the airplanes is not 100% and this is accounted for by correcting wjk to w ˜ jk = wjk (α + η(1 − α)) where α = 0.7 corresponds to the average occupancy rate of order of 70% provided by IATA and η is a random number drawn uniformly in the interval [−1, 1] at each time step

Influenza – Colizza et al. (2007b) Cities

Influenza – Colizza et al. (2007b) Countries

Influenza – Colizza et al. (2007b) Continents

Influenza – Colizza et al. (2007b)

Movies

Individual-based models (IBMs) IBMs treat individuals as separate modeling units, and follow their fate from their birth to their death SIS model with no demography: consider a lattice composed of square sites, each of which is occupied by either a susceptible or an infectious individual and whose potential transitions are governed by their neighbors: S →I

at rate β × number of I in the neighborhood of S

I →S

at rate γ

Lattice initially full of susceptibles, only some individuals are inoculated by a parasite

Individual-based models (IBMs)

β = 2, γ = 1

Individual-based models (IBMs) Spatial IBMs are simulation models, involving some stochasticity and producing plenty of numerical output – rigorous analysis far from feasible once they become a bit more complex Simple IBMs are an abstraction of reality but superb tools to understand how the individual and spatial nature of populations causes epidemic dynamics to deviate from their (deterministic and) globally interacting counterparts IBMs can also be formulated in discrete time and/or continuous space and can be designed to include a variety of complex and detailed host behavior This can however lead to a huge number of parameters or functional relationships that can be difficult to determine from available data

Network approaches It is contact structure between individuals that determines disease progression through the host population One of the simplest assumptions is that contact structure forms a network of links between individuals (or nodes), with all links of equal strength

Network approaches

Lattice-like IBMs: rigid spatial structure where each individual has the same number of local neighbors Network models are far more flexible and allow individuals to have different numbers of contacts some of which may be distant For airborne diseases, such as influenza or measles, deciding who is a (social) contact is very difficult though still important For sexually transmitted diseases, it is far clearer how to define (sexual) contacts even if data is somewhat more delicate

Network approaches

The primary advantage of network models is their ability to capture complex individual-level structure in a simple framework All interaction strengths are usually considered identical, interactions symmetric and networks undirected (infection can pass in both directions across a contact, with equal probability) Types of commonly used networks differ by such fundamental features as the heterogeneity in the number of contacts, clustering of contacts, and average path length (number of steps it takes on average to link two randomly chosen individuals)

Network approaches

Random network – probability of any two individuals being connected is the same Spatial network – probability of any two individuals being connected depends on the distance between them Scale-free network – most individuals have relatively small number of contacts and a few have many contacts (each new individual that is added to the population connects preferentially with individuals that already have a larger number of contacts)

Network approaches

random, spatial and scale-free networks

Network approaches

Small-world network – basically a lattice structure, with a small number of long-range connections added

lattice and small-world networks

Network approaches

Networks have many similarities to lattice-based IBMs, so the simulations proceed along the same lines SIS model with no demography S →I

at rate β × number of infectious neighbors

I →S

at rate γ

Network approaches

random, spatial and scale-free networks (β = 1, γ = 0.5)

Network approaches

lattice and small-world networks (β = 1, γ = 0.5) Solid to dashed lines in the right panel: increasing the number of long-range contacts from 10 to 20 and 100

Network approaches Network models display slower epidemic dynamics compared with random-mixing models Network models that are most like the random-mixing models – with short average path length (small-world, random and scale-free networks) and little clustering (random and scale-free networks) – demonstrate the fastest epidemic dynamics for a given average number of contacts per individual Network models have applications to many communicable diseases and are particularly appropriate for sexually transmitted diseases where the network structure may be far more apparent Much remains to be done, both on the theoretical side (approximations, network models are computationally quite intensive) and the applied side (parameterization)

Summary: which spatial model to use?

All spatial models share a common theme: local interactions generally dominate longer-range interactions, leading to localized disease transmission and hence clustering of cases Selected model type should reflect the host species, the problem being addressed (including spatial scale of interest), availability of data, and the form of required results (detail and accuracy of predictions) Model choice is a skill in itself: ideally, several model types should be developed, their output tested to the available data, and their predictions scrutinized in terms of questions addressed and detail of answers required

Infectious diseases and control

Infectious diseases and control

Infections is something we virtually always perceive as undesirable and harmful We would like to avoid any disease epidemic, and if an infection is already present, we would like to minimize its impact and eliminate it altogether, at best Approaches applied to avoid or mitigate impacts of an infection include vaccination, treatment, isolation and quarantine, age-structured culling, prevention, and vector control Infections can also be welcome: control of plant and animal pests

Control of infectious diseases

Threshold behavior of many epidemiological models: R0 < 1 = infection will die out R0 > 1 = there will be an epidemic R0 might be a complex function of model parameters An obvious way to prevent epidemics is to get and keep R0 below 1, in one way or another When an infection is already present in the population, the same holds for the (effective) reproduction number R

Control of infectious diseases

Epidemic SIR model Mass action incidence (βSI ): R0 = βN/γ ⇒ there is a host density threshold NT below which the pathogen cannot invade population of (susceptible) hosts Standard incidence (βSI /N): R0 = β/γ ⇒ no such threshold Culling to reduce the (susceptible) host population below NT is a common policy for handling disease outbreaks in wildlife or in domestic animals (fails if there is no host-density threshold)

Control of infectious diseases

Endemic models = stable positive equilibrium = steady disease prevalence This equilibrium is a complex function of model parameters so that manipulating any of them in a proper way can result in a decrease in disease prevalence

Vaccination

Herd immunity = enough people have disease-acquired or vaccination-acquired immunity, so that the introduction of one infective into the population does not cause an invasion of the disease The ultimate goal of any vaccination program = vaccinate enough people to achieve herd immunity Successful vaccination moves individuals from the susceptible class (newborns or older ones) straight into the recovered class so that they can no longer catch or spread the infection

Vaccination

Vaccination of newborns (p is the fraction of newborns vaccinated): SI dS = bN(1 − p) − β − dS dt N SI dI = β − γI − dI dt N dR = bNp + γI − dR dt Analysis (b = d): disease cannot invade if p > pc = 1 − (γ + d)/β = 1 − 1/R0

Vaccination

Disease Smallpox Measles (rural USA) Measles (rural USA) Measles (UK)

R0 5 5.4-6.3 8.3-13.0 12.5-16.3

pc [%] 80 81.5-84.1 88.0-92.3 92.0-94.0

Smallpox: the only disease for which this has actually been achieved worldwide and smallpox has been eliminated (the last known case was in Somalia in 1977, and the virus is maintained now only in laboratories) Measles: vaccine is not always effective, and vaccination campaigns are never able to reach everyone → herd immunity against measles has not been achieved (and probably will never be)

Vaccination

Vaccination Vaccination of susceptibles (θ is the rate at which susceptibles are vaccinated; vaccination of I and R individuals is assumed to have no effect): SI dS = bN − β − dS − θS dt N SI dI = β − γI − dI dt N dR = θS + γI − dR dt Analysis (b = d): disease cannot invade if θ > θc = d(β/(γ + d) − 1) = d(R0 − 1)

Treatment Treatment for infection Fraction γ per unit time of infectives are selected for treatment Treatment reduces infectivity by fraction δ Epidemic SITR model (T is the treatment class): dS/dt = −βSI /N − δβST /N

dI /dt = βSI /N + δβST /N − γI − αI dT /dt = γI − ηT

dR/dt = αI + ηT

I (t) → 0 and T (t) → 0 as t → ∞

Treatment Basic reproduction number γ δβ β + α+γ α+γ η

2

2

1.5

1.5

RT

RT

RT =

1

0.5

0 0

1

0.5

0.5

1 η

1.5

2

0 0

0.5

1 δ

1.5

2

Quarantine and isolation For an outbreak of a new disease such as SARS, to which no vaccine or treatment are available, quarantine and isolation are the only control measures available dQ . . . rate at which infectives are detected and quarantined 1/τ . . . average time spent in quarantine (individuals are assumed to leave the quarantine class only after they have recovered) Epidemic SIQR model (Q is the quarantine class): dS/dt = −βSI /N

dI /dt = βSI /N − γI − dQ I

dQ/dt = dQ I − τ Q

dR/dt = γI + τ Q

Quarantine and isolation

Basic reproduction number RQ =

β γ + dQ

RQ < 1 if and only if dQ > dQ∗ = β − γ More complex model: exposed members (imperfectly) quarantined and diagnosed infectives and quarantined members (imperfectly) isolated

Vector control For infections transmitted by vectors, measures to control vector populations essentially stand for measures to directly control the infection This has been attempted, e.g. for mosquitoes by inundative releases of sterile mosquito males to limit mosquito reproductive potential and hence to mitigate impacts of malaria and other mosquito-borne infections (yellow fever, dengue, . . .) Vector-host models: R0 < 1 if and only if N βs βm > Z b+γ c If the ratio N/Z of the total host population to the total vector population becomes large, the infection goes extinct

Control by infectious diseases

Pathogens also considered as control agents of some pest species Question: Can a pathogen exert long-term control of its host, and more specifically can it regulate the abundance of a pest below its economic threshold? Ethical and safety questions are another side of the story Pathogens may cause an extra mortality or reduce natality of infected hosts We might expect a decline in the host population growth rate or a decrease in a stable endemic equilibrium of the host

Control by infectious diseases Exponentially growing host populations Parasites can reduce the growth rate or even regulate such populations to a stable (including zero) equilibrium 40 35

40 no disease disease, α > α

35

disease, α < αT

30

30

Total host population

Total host population

T

25 20 15 10 5 0 0

ˆ 1, R1 > 1 disease, R ˆ > 1, R1 < 1 disease, R

25 20 15 10 5

50

100

150 Time

200

250

300

0 0

20

40

60 Time

80

100

Control by infectious diseases Logistically growing host populations Looking at how the endemic equilibrium depends on disease-related parameters, we can in principle suggest pathogen properties that would be optimal for achieving pest control goals

12 Total host population

Total host population

12 10 8 6 4 no disease; disease, R0 < 1

2 0 0

disease, R0 > 1 200

400

600 Time

800

1000

10 8 6

no disease; disease, R < 1

4

disease, R0 > 1, i* < (b−d)/α

2

disease, R0 > 1, i* > (b−d)/α

0 0

0

500

1000 Time

1500

2000

VVIC (virus-vectored immunocontraception)

Pest populations inoculated by a virus that is able to sterilize either males or females, but is able to spread through both sexes These VVIC agents are often genetically modified forms of natural viral forms (such as myxomatosis in rabbits) Assumptions: promiscuous mating system, common to mammals that are the prime target for VVIC; infection directly transmitted between animals by sufficiently close contacts

VVIC (virus-vectored immunocontraception)

Pre-control phase each male is assumed to mate with h randomly chosen females within each reproductive cycle each female that mates with at least one male gives birth to b offspring, each becoming male or female with probability 0.5 only a fraction 1 − N/K of the offspring survive to adulthood (N is the current adult population size and K is the carrying capacity of the environment) each adult individual dies with probability d at the end of reproductive cycle

VVIC (virus-vectored immunocontraception)

Control phase Non-mating period infection is transmitted between random pairs of individuals through Λ(N) non-mating encounters if only one of the pair members bears the virus, the other becomes infected with probability β Mating period each male mates with h randomly chosen females upon mating, infection is also transmitted with probability β

VVIC (virus-vectored immunocontraception)

Control phase Reproduction male-specific virus: any infected male has probability 1 − σ to successfully fertilize a female upon mating, and one successful mating is assumed sufficient for the female to produce b offspring female-specific virus: mated females produce b offspring if susceptible, and b or no offspring with probability 1 − σ and σ, respectively, if infected Mortality as in the pre-control phase

VVIC (virus-vectored immunocontraception)

VVIC (virus-vectored immunocontraception)

Viruses used as VVIC agents are often genetically modified forms of natural viral forms (such as myxomatosis in rabbits) These two forms compete for susceptible hosts This requires consideration of multi-strain interactions, yet single-strain models are useful as an a priori assessment of the VVIC efficiency provided the GM virus is competitively superior

Control by infectious diseases

Different perspective on control by infectious diseases Rapid population growth in the eighteenth century: – Europe: from 118 million to 187 million – Great Britain: from 5.8 million to 9.15 million – China: from 150 million to 313 million Plausible explaination: – improvements in agriculture and food production – decrease in the death rate due to diseases

Control by infectious diseases

Disease death rates dropped sharply in the eighteenth century, partly from better understanding of the links between illness and sanitation and partly because the recurring invasions of bubonic plague subsided, perhaps due to reduced susceptibility. One plausible explanation for these population increases is that the bubonic plague invasions served to control the population size, and when this control was removed the population size increased rapidly.

Disease evolution

Disease evolution

Parasite action on the host will select for individuals with reduced susceptibility to the infection Selective forces also act strongly on parasites – persistence of an infection is facilitated by low pathogenicity and by long duration of the infection Pathogens may constitute the selective force responsible for the evolution and maintenance of sexual reproduction in animal and plant species – the Red Queen hypothesis In natural populations, variability in host susceptibility to an infection appears to be the rule rather than the exception

Disease evolution

Types of evolution commonly studied Evolution of (disease) virulence Evolution of host resistance Co-evolution of (disease) virulence and host resistance Evolution of drug-resistant viral strains

Disease evolution The myxoma virus and rabbits in Australia

Disease evolution

The myxoma virus and rabbits in Australia

Evolution of virulence

Classic approach: maximization of R0 as a function of virulence Population-level perspective Rationale: Since R0 is a measure of disease transmissibility, maximizing a pathogen’s R0 is equivalent to maximizing its transmissibility and (initial) spread rate through the host population

Alternative approach: invasion analysis Individual-level perspective

Evolution of virulence

Two virus strains R0 maximization: both start in the pristine environment – who spreads faster? Invasion analysis: one strain is resident, the other strain invades an established resident strain-host endemic equilibrium – does the invader succeed and expel the resident? R0 maximization works for relatively simple epidemiological models, but it doesn’t work in general R0 maximization is relatively easy so it is good to know when the two approaches produce identical results and when they deviate

Evolution of virulence

SIS model x . . . disease virulence due to a virus strain dS = bS (x, S, I )S + bI (x, S, I )I − β(x, S, I )SI − d(x, S, I )S + γ(x, S, I )I dt dI = β(x, S, I )SI − d(x, S, I )I − α(x, S, I )I − γ(x, S, I )I dt

Evolution of virulence

Invasion analysis Let the resident population with strain virulence x be at endemic equilibrium (S ∗ (x), I ∗ (x)) Mutant strain with virulence x 0 is now introduced into the resident population Mutant strain is initially rare I 0 . . . density of hosts infected by strain x 0 Assumption: no super-infection or co-infection

Evolution of virulence

Growth rate of I 0 class: dI 0 = f (x, x 0 )I 0 dt f (x 0 , x) = mutant’s per capita growth rate or invasion fitness f (x, x 0 ) = β(x 0 , S ∗ (x), I ∗ (x))S ∗ (x) − d(x 0 , S ∗ (x), I ∗ (x)) −α(x 0 , S ∗ (x), I ∗ (x)) − γ(x 0 , S ∗ (x), I ∗ (x))

Evolution of virulence

Evolutionarily stable equilibrium (ESS) is attained for a virus strain x ∗ for which the selection gradient g (x) ≡

∂ f (x 0 , x)|x 0 =x = 0 ∂x 0

No mutant strain with lower or higher virulence than x ∗ can invade the resident strain with virulence x ∗ The resident strain with virulence x ∗ cannot be invaded by any mutant strain

Evolution of virulence

Reproduction number of mutant strain with virulence x 0 R(x 0 , x) =

β(x 0 , S ∗ (x), I ∗ (x))S ∗ (x) d(x 0 , S ∗ (x), I ∗ (x)) + α(x 0 , S ∗ (x), I ∗ (x)) + γ(x 0 , S ∗ (x), I ∗ (x))

Mutant can invade if f (x 0 , x) > 0 ⇔ R(x 0 , x) > 1

Evolution of virulence R(x, x) = 1 . . . density of infectives accurately replenishes itself once the disease with resident stain x has reached an endemic equilibrium This implies S ∗ (x) =

d(x, S ∗ (x), I ∗ (x)) + α(x, S ∗ (x), I ∗ (x)) + γ(x, S ∗ (x), I ∗ (x)) β(x, S ∗ (x), I ∗ (x))

Hence R(x 0 , x) =

β(x 0 ,S ∗ (x),I ∗ (x)) d(x 0 ,S ∗ (x),I ∗ (x))+α(x 0 ,S ∗ (x),I ∗ (x))+γ(x 0 ,S ∗ (x),I ∗ (x)) β(x,S ∗ (x),I ∗ (x)) d(x,S ∗ (x),I ∗ (x))+α(x,S ∗ (x),I ∗ (x))+γ(x,S ∗ (x),I ∗ (x))

Evolution of virulence Basic reproduction number of resident strain R0 (x) =

β(x, S0 , 0)S0 d(x, S0 , 0) + α(x, S0 , 0) + γ(x, S0 , 0)

Basic reproduction number of mutant strain R0 (x 0 ) =

β(x 0 , S0 , 0)S0 d(x 0 , S0 , 0) + α(x 0 , S0 , 0) + γ(x 0 , S0 , 0)

When parameters d, β, α and γ are density-independent: R(x 0 , x) = R0 (x 0 )/R0 (x) Therefore R(x 0 , x) > 1 ⇔ R0 (x 0 ) > R0 (x) which implies that invasion analysis and R0 maximization produce equal results

Evolution of virulence

Example 1: Evolution towards no harm Assumptions: bS (x, S, I ) = bI (x, S, I ) = b, β(x, S, I ) = β, d(x, S, I ) = d, γ(x, S, I ) = γ, b = d Let α(x, S, I ) = x . . . virulence = disease-induced host mortality S ∗ (x) = (x + d + γ)/β and hence f (x, x 0 ) = x − x 0 f (x 0 , x) > 0 ⇔ x 0 < x ⇒ mutant strain with lower virulence than resident strain can always invade and therefore the system will evolve to the most benign strain

Evolution of virulence Example 1: Evolution towards no harm The same conclusion is achieved by maximizing R0 (x) = βS0 /(x + d + γ) with respect to x 6 5

R0

4 3 2 1 0 0

0.5

1

1.5 x

2

2.5

3

Evolution of virulence The virulence-transmission trade-off Disease-related parameters are likely not uncorrelated Rather, they are often likely to be traded off Example: disease-induced mortality α and transmission efficiency β are often assumed related through a hyperbolic function such as β(α) = α/(α + c), for some c > 0 1 0.8

β

0.6 0.4 0.2 0 0

2

4

α

6

8

10

Evolution of virulence Trade-off between transmissibility and disease-induced mortality An increase in infectiousness β can only be achieved at the expense of a shortening of the longevity of infectious hosts Think of α and β as functions of the number of copies of the parasite circulating in an infected individual (parasitemia): The number of copies is small: the chance of new infection is small, and mortality from the disease of an infected individual is likewise small The number of copies raises: the chance of transmission also raises, but eventually levels off because, for example, the parasite is saturating the exhaled air of an infected individual; mortality from the disease of an infected individual raises as well

Evolution of virulence

Example 2: Evolution to an intermediate virulence Assumptions as above Invasion fitness f (x 0 , x) = (x − x 0 )[x x 0 − c(d + γ)]/[x(x 0 + c)] ⇒ evolution attains an intermediate virulence x ∗ =

p

c(d + γ)

Evolution of virulence Example 2: Evolution to an intermediate virulence The same conclusion is achieved by maximizing R0 (x) = x/(x + c) S0 /(x + d + γ) with respect to x 2

R0

1.5

1

0.5 *

x 0 0

0.5

1

1.5 x

2

2.5

3

Evolution of virulence

Example 3: Evolution need not maximize R0 γ = γ0 /(1 + I /K ) for some K > 0 . . . care imposed to infected hosts declines with their overall density Other assumptions as in example 2 b = 2, d = 1, c = 1, γ0 = 1, K = 10 Invasion analysis: x ∗ = 1.061 (also x ∗ = 1.253 for b = 1.75 and x ∗ = 1.367 for K = 100)

Evolution of virulence Example 3: Evolution need not maximize R0 R0 (x)p = x/(x + c) S0 /(x + d + γ0 ) is maximized at x˜∗ = c(d + γ0 ) ⇒ no effect of K and b √ ⇒ x˜∗ = 2 ≈ 1.414 2

R0

1.5

1

0.5 x* 0 0

0.5

1

\tilde{x}^* 1.5 x

2

2.5

3

Both quantitative and qualitative effects missed when maximizing R0

¡odels of macroparasitic infections

Microparasites versus macroparasites

Microparasites viruses, bacteria, (parasitic) protozoa (or protozoans), and fungi Macroparasites (parasitic) helminths, some arthropods such as ticks and mites, and (parasitic) protozoa (or protozoans)

Microparasites versus macroparasites

Microparasites small size – generally cannot be seen with the naked eye Macroparasites large enough to be seen with the naked eye

Microparasites versus macroparasites

Microparasites duration of infection typically short, of a transient nature Macroparasites infections tend to be of persistent nature, with hosts being continually reinfected

Microparasites versus macroparasites

Microparasites tendency to induce immunity to reinfection in those hosts that survive the initial onslaught Macroparasites immune responses elicited generally depend on the number of parasites present in a given host, and tend to be of a relatively short duration

Microparasites versus macroparasites

Microparasites short generation times Macroparasites much longer generation times than microparasites

Microparasites versus macroparasites

Microparasites extremely high rates of reproduction within the host Macroparasites multiplication within the host either absent or occurs at low rates

Microparasites versus macroparasites

Microparasites can complete full life cycle inside a single host Macroparasites grow in the host but multiply by producing infective stages such as eggs, spores or cycts that, as a developmental necessity, are released from the host

Differences in the life cycles of microparasites versus macroparasites and impacts on their hosts require development of different epidemiological models

Models of macroparasitic infections Birth and/or death rates of infected hosts, parasite virulence, rate of production of transmission stages and any host resistance to infection all typically depend on numbers of parasites these hosts harbor Variability in worm burden → one must track numbers of parasites (or parasite burden) within each host Modelers distinguish categories of hosts, with p(i ) denoting the frequency of hosts with i parasites, i = 0, 1, 2, . . . (p(0) is the frequency of healthy hosts) Hi = p(i )H is then the number of hosts harboring i parasites (H = total host population size)

Models of macroparasite infections

Basic model

H(t) = host population size at time t P(t) = (within-host) parasite population size at time t P/H = average number of parasites per host p(i ) = frequency of hosts with i parasites = probability that a host contains i parasites → will depend on the probability distribution of parasites within the host population

Basic model

Host population growth rate (in the absence of infection) – no density-dependent constraints on population growth assumed bH − dH Assumption: b > d Without parasites, the host population grows exponentially

Basic model Parasite-induced host mortality rate – assumed linearly proportional to the number of parasites a host harbors α . . . virulence of one parasite in a host αi . . . parasite-induced mortality rate of a host with i parasites Total rate of loss of hosts ∞ ∞ X X ip(i ) αiHi = αH i =0

i =0

Because ∞ X

ip(i ) = mean parasite burden = E [i ] = P/H

i =0

⇒ parasite-induced host mortality rate is αP

Basic model

Parasite fecundity and transmission λ . . . rate of production of transmission stages per parasite Rate of production for the total parasite population ∞ X i =0

λiHi = λH

∞ X i =0

ip(i ) = λP

Basic model

Parasite fecundity and transmission Direct life cycle (only a single host species involved): transmission stages will pass out of the host into the external habitat and will survive there as resistant stages of free-living larvae awaiting contact with or ingestion by another host In the external habitat: transmission stages are subject to natural mortality due to senescence or predation → only a proportion of those released successfully enters new hosts → this proportion depends on the size of the host population relative to other “absorbers” of transmission stages

Basic model Parasite fecundity and transmission Assumption: transmission virtually instantaneous = no time delay due to transmission or development between the birth of a transmission stage and successful contact with a new host (this holds for some parasites but in the majority certain developmental processes have to occur before the transmission stage becomes fully infective) Simple model: transmitted proportion = H/(H + H0 ) (H0 constant) Rate at which new parasites are acquired by hosts λP

H H + H0

Basic model Parasite mortality rate – three components: 1. Losses due to natural host mortality d: ∞ X i =0

dHi i = dH

∞ X

ip(i ) = dP

i =0

2. Losses from parasite-induced host deaths: ∞ ∞ X X i 2 p(i ) = αHE [i 2 ] (αi )Hi i = αH i =0

i =0

3. Losses due to natural parasite mortality µ within the host (deaths due to host immunological responses and parasite senescence): µP

Basic model Dynamics Rate of change of the host population: dH = (b − d)H − αP dt Rate of change of the parasite population: H dP = λP − dP − µP − αHE [i 2 ] dt H + H0 Problem: E [i 2 ] – precise value depends on the distribution p(i ) Simplification: a phenomenological assumption about the statistical distribution of parasites among hosts

Basic model

Aggregated (or over-dispersed): negative binomial parasite distribution p(i ) =



k  i   k +i −1 k k 1− k + P/H k + P/H i

Parameter k provides an inverse measure of the degree of parasite aggregation within the host population – the smaller k the greater the degree of parasite clumping 2

E [i ] =



P H

2

k +1 P + k H

Basic model

(Independently) random: Poisson parasite distribution – parasites distributed independently randomly within hosts (corresponds to k → ∞ in the negative binomial distribution) p(i ) =

2

(P/H)i −P/H e i!

E [i ] =



P H

2

+

P H

Basic model Regular (or under-dispersed): positive binomial parasite distribution Parasites are distributed more regularly than randomly, which technically corresponds to replacing k by −k 0 in the negative binomial distribution; parameter k 0 is now inversely proportional to the degree of regularity – in the limit k 0 → 0, all hosts carry identical parasite burden

p(i ) =



−k 0 + i − 1 i

 

2

E [i ] =

−k 0 −k 0 + P/H



P H

2

−k 0 

k0 − 1 P + k0 H

−k 0 1− −k 0 + P/H

i

Basic model Almost without exception observed patterns are aggregated, where a relatively few members of the host population harbor the majority of the total parasite population → negative binomial distribution has proved to be a good empirical model for a large number of observed parasite distributions

Basic model

Basic model Poisson distribution of p(i ): E [i 2 ] = (P/H)2 + P/H dH = (b − d)H − αP dt dP H P2 = λP − (d + µ + α)P − α dt H + H0 H Equilibrium P ∗ /H ∗ = (b − d)/α H ∗ = H0

µ+α+b λ − (µ + α + b)

Parasites are able to regulate hosts if λ>µ+α+b

Basic model Stability analysis: equilibrium (H ∗ , P ∗ ) is neutrally stable

Basic model

λ < µ + α + b → host (and also parasite) population will grow exponentially; mean number of parasites per host will tend to zero Basic model is structurally unstable: the slightest alteration in the form of model assumptions will precipitate the system either to stability (disturbances damping back to the equilibrium) or to instability (growing oscillations out of the equilibrium) Structural instability makes the model biologically unrealistic! But: Basic model is a very useful tool as a point of departure for various biologically realistic modifications, as it helps reveal mechanisms that stabilize or destabilize host-parasite dynamics

Negative binomial distribution

E [i 2 ] = Model

 P 2 k+1 H k

+

P H

dH = (b − d)H − αP dt H P2 k + 1 dP = λP − (d + µ + α)P − α dt H + H0 H k

Equilibrium P ∗ /H ∗ = (b − d)/α

H ∗ = H0

µ + α + d + (b − d)(k + 1)/k λ − (µ + α + b + (b − d)(k + 1)/k)

Negative binomial distribution

Parasites are able to regulate hosts if λ > µ + α + d + (b − d)(k + 1)/k Stability analysis: equilibrium (H ∗ , P ∗ ) is globally stable (damping oscillations) for any finite k > 0 Reverse condition: the host population grows exponentially (as above) Positive binomial distribution: no stable equilibrium

Negative binomial distribution

Negative binomial distribution

Negative binomial distribution

Many parasites highly aggregated (small k) For the parasite to control the host population, α needs to be relatively low and/or λ relatively large (nematode Haemonchus contortus, a parasite of sheep – 10000 eggs per day) Exact form of the aggregated distribution is immaterial – negative binomial distribution can be replaced with other qualitatively similar distributions with the same conclusions

Non-linear parasite induced mortality rate

Non-linear parasite induced mortality rate Assumption: αi → αi 2 = parasite-induced mortality rate of a host with i parasites Model

dH = (b − d)H − αHE [i 2 ] dt H dP = λP − (d + µ)P − αHE [i 3 ] dt H + H0

Enhancement in equilibrium stability, relative to the previous models: → (H ∗ , P ∗ ) stable for aggregated and Poisson distributions → (H ∗ , P ∗ ) stable for regular dist. if k 0 is not too small

However: → domain of parameter space for which (H ∗ , P ∗ ) exists is smaller

Non-linear parasite induced mortality rate

In words: → more steeply severe parasite-induced host mortality makes it relatively harder for the parasite to check the host population growth → but if it can do so, this kind of density dependence results in relatively faster damping of perturbations and enhanced equilibrium stability Non-linear parasite induced mortality rate helps stabilize the host-parasite system

Density dependence in parasite population growth

Severity of parasite-induced immunological responses or intraspecific competition for finite resources such as space or nutrients increases with within-host parasite burden → parasite mortality rate µ increases with parasite burden i and/or

→ parasite reproductive rate λ decreases with parasite burden i

Density dependence in parasite population growth

Density dependence in parasite population growth Assumption: µ → µi 2 = natural parasite mortality within the host Model

dH = (b − d)H − αP dt H dP = λP − dP − (µ + α)HE [i 2 ] dt H + H0

Enhancement in equilibrium stability, relative to the previous models: → (H ∗ , P ∗ ) stable for aggregated and Poisson distributions → (H ∗ , P ∗ ) stable for regular dist. if k 0 is not too small

However: → domain of parameter space for which (H ∗ , P ∗ ) exists is smaller

Density dependence in parasite population growth λ a decreasing function of i

Modeling: qualitatively similar dynamical properties to those for density dependent mortality rate µ Density dependence in parasite growth helps stabilize the host-parasite system

Parasite-induced reduction in host reproduction Effects commonly produced by helminths or arthropod parasites, in particular by larval helminths within invertebrate intermediate hosts

Assumption: b → b − βi = per capita reproductive rate of hosts P Averaging over hosts: b → b − β ∞ i =0 ip(i ) = b − βP/H

Parasite-induced reduction in host reproduction Crude model: host’s reproductive rate cannot be negative But: captures the essentials of the phenomenon (and enables mathematically elegant calculations) dH = (b − d)H − (α + β)P dt H dP = λP − (µ + d)P − αHE [i 2 ] dt H + H0 Aggregated distribution: regulation if λ > µ + α + d + [α(b − d)(k + 1)]/[k(α + β)] equilibrium globally stable if α > βk Random and regular dist.: no stable equilibrium possible

Parasite-induced reduction in host reproduction

Parasite-induced reduction in host reproduction

Stability seriously impaired if the parasite has a significant influence on the host’s reproductive rate, and possible only within a narrowly constraint domain of parameter space (β and/or k small, but small k makes difficult to satisfy the other condition) Many natural host-parasite associations persist, despite the parasite substantially altering the host’s reproductive rate – parameters characterizing those natural associations are likely a far-from-random set, or other, stabilizing mechanisms are also present

Parasites reproducing within their hosts

Many parasites such as various protozoans reproduce within their final host, directly contributing to the growth of their population within their host Distinct process from production of transmission stages Two types of reproduction thus occur r . . . rate (per parasite) at which parasites reproduce within the host

Parasites reproducing within their hosts Model

dH = (b − d)H − αP dt H dP = λP + rP − (d + µ)P − αHE [i 2 ] dt H + H0

Aggregated distribution: (H ∗ , P ∗ ) globally stable if λ + r > µ + α + d + (b − d)(k + 1)/k > r r > µ + α + d + (b − d)(k + 1)/k . . . direct parasite reproduction is so great as to produce excessively severe host mortality → both species go extinct = new type of dynamical behavior λ + r < µ + α + d + (b − d)(k + 1)/k . . . parasite population cannot increase fast enough to control hosts → both populations grow exponentially

Parasites reproducing within their hosts

Parasites reproducing within their hosts

Large r (real r tend to be large): extinction may be avoided if k is small Random distribution: equilibrium is (globally) neutrally stable Regular distribution: equilibrium (if it exists) is always unstable Destabilizing mechanism. The fact that many natural host-parasite associations persist, despite the parasites’ having a high rate of direct reproduction suggests that other interaction parameters have coevolved in such a way as to maintain the continued existence of both populations, or other, stabilizing mechanisms are also present

Effect of time delays

All the preceding models: no time delays assumed to occur between production of a transmission stage and its availability for reinfection of a new host – eggs, spored or cysts produced by the parasite assumed sufficiently developed to enable them successfully infect a host upon a contact Many parasites exhibit developmental delays between the point of production of a transmission stage and its readiness for reinfection (actually, the length of such delays will often be associated with climatic factors influencing the habitat) Third ODE needed that gives an explicit account of population dynamics of the free-living infective stage W

Parasites reproducing within their hosts

Effect of time delays Assumption: aggregation of parasites within hosts dH = (b − d)H − αP dt dP P2 k + 1 = βWH − (µ + α + d)P − α dt H k dW = λP − γW − βWH dt β . . . rate at which hosts pick up infective stages (∼ mass action incidence) γ . . . rate of loss of infective stages from the habitat (due to death or other processes that foreclose possibility of their infecting a host)

Effect of time delays

(H ∗ , P ∗ , W ∗ ) exists if λ > µ + α + d + (b − d)(k + 1)/k

(H ∗ , P ∗ , W ∗ ) is stable if (QUITE COMPICATED CONDITION) If λ, γ and β (birth, death and transmission rates of W ) are high relative to other rates (often the case), approximate stability condition is 1 µ + α + d + (b − d)(k + 1)/k . γ k . . . otherwise populations oscillate in a stable limit cycle (hence equilibrium is unstable)

Effect of time delays

As 1/γ (mean time of persistence of infective stages in the habitat if there is no host picking) becomes large, more extreme aggregation (smaller k) is required for the host population to be regulated by parasites (that is, for the equilibrium to be stable) Time delays in the development or transmission of infective stages of the parasite have a destabilizing influence (condition for equilibrium existence same as without W , but for its stability much more rigid)

Summary

Examined modifications of the basic (neutrally stable) model Stabilizing forces (Anderson and May 1978) aggregated parasite distribution among hosts non-linear parasite-induced mortality rate of hosts various kinds of density dependence in parasite growth Destabilizing forces (May and Anderson 1978) parasite-induced reduction in host reproduction parasites directly reproducing within hosts time delays in parasite development and transmission

Summary

Real host-parasite associations likely exhibit all the above effects, to a greate or lesser extent – there is a tension between the stabilizing and destabilizing elements of the interaction Modeling shows that parasite species are capable of regulating host population growth, even in the absence of other influences such as predation or intraspecific competition among hosts To examine the (real) importance of parasites in affecting the growth characteristics of their host population, observations are ideally required for populations of hosts with and without parasite infection

Summary

Some extensions host population assumed constant or growing in a density-dependent way hypobiosis – consumed infective stages of the parasite are dormant in the host for some time before becoming adults (Dobson and Hudson 1992, nematode Trichostrongylus tenuis in red grouse) macroparasites with indirect life cycles (primary and intermediate hosts) – schistosomiasis (man and snail) sexual reproduction of some macroparasites (majority of worms reproduce sexually within humans → mating function → possibility of a mate-finding Allee effect)