We now have the basic framework and toolkit needed to start thinking about building an infectious disease model. As you may have already realised, much of the hard work lies in conducting a thorough literature review to develop a sound understanding of how the disease operates in the target population. The next step in translating veterinary knowledge into a simulation model is to develop a conceptual representation of how animals become infected and what happens to them as they progress through different disease states. This includes identifying relevant health and production impacts that may be tracked to evaluate the costs of disease and the potential benefits of control. In this module, we focus on developing simple discrete-time compartmental models that form the foundation for more sophisticated within-herd and between-herd modelling approaches.
By the end of this module, students should be able to:
There are many different ways we can choose to represent the transmission and pathogenesis of infectious disease. Compartmental models are one of the oldest and most widely used approaches to modelling infectious disease dynamics in populations with the simplest being the so-called SIR model.
In this model, all individuals in the population of size N are assumed to belong to one of the three mutually exclusive disease states. Disease states are represented by boxes (compartments) and the progression of individuals between the compartments is shown with arrows. In this case, Susceptible (S) individuals become Infectious (I) at a rate proportional to both the number of infectious individuals in the population and the transmission coefficient Beta (β), which is a combined estimate of the average contact rate between individuals and the probability of disease transmission occurring through the contact. You can think of βI as representing the infection pressure on individuals in the population, which we would expect to increase as the number of infected individuals in the population also increases. Once infected, individuals become Recovered (R) and immune to the disease at the rate Gamma (γ), which is estimated as the inverse of the infectious period (i.e. if the average infectious period lasts 7 days, we would expect 1 individual to move out of the I compartment every 7 days on average = 1 individual / 7 days or 0.14 individuals per day). We will talk about methods for estimating β in Module 6 where we go into more detail on infectious disease dynamics.
In an epidemic situation with a population size N, we would typically start with 1 individual in the I compartment as the index case, (N-I) individuals in the S compartment as an entirely susceptible population, and 0 individuals in the recovered compartment. For diseases that operate over longer time scales, we can also build in parameters to describe population turnover through birth rates Mu (μ) that replenish the population of S individuals and death rates Alpha (α) that describes the death rates from natural processes. In a steady state system with no change in population size, the birth rate and death rate will be equal and can therefore be represented as a single population turnover rate (μ).
There are almost infinite ways we can modify the basic SIR model to reflect our disease system of interest. We typically use capital English letters to name each compartment and Greek symbols to represent the different transition rates between compartments. Although there are some existing conventions around which letters and symbols are used to represent common process (i.e. using β to describe transmission and E to describe individuals who have been exposed to the disease, but are not yet infectious), to a large extent it doesn’t really matter what you call them as long as everything is defined clearly in your write-up. This is best done by providing readers with a diagram similar to the above figures along with a summary table similar to the one below. For each parameter, it is important to provide a clear definition of what it represents as well as the baseline values you are using in the model, the units for the values, and a reference for where you obtained the estimate for the values. Your text should also include a written summary of the disease dynamics to help readers understand the rationale behind the model.
| Parameter | Definition | Baseline value | Units | Reference |
|---|---|---|---|---|
S |
Susceptible individuals | 95 | Number of individuals | - |
I |
Infectious individuals | 5 | Number of individuals | - |
R |
Recovered individuals | 0 | Number of individuals | - |
N |
Total population size | 100 | Number of individuals | - |
β |
Transmission coefficient | 0.007 | Rate per day | Smith et al |
γ |
Recovery rate | 0.14 | Individuals per day | Smith et al |
μ |
Birth rate | 0.01 | Individuals per day | Expert opinion |
α |
Death rate | 0.01 | Individuals per day | Expert opinion |
By now, you might be wondering how we actually use these compartmental models to start understanding the disease transmission dynamics. Initially, we want to be able to track the number of individuals in each compartment over time so that we can see how fast disease is spreading through the population. To do this, we need to turn the conceptual model into a series of simple differential equations to describe the change in the number of individuals in each compartment per unit time.
The number of equations will always be equal to the number of compartments and the number of terms within each equation will always be equal to the number of arrows entering or exiting the compartments. Terms describing arrows in will be positive since they are adding individuals to the compartment while terms describing arrows out will always be negative since they are removing individuals from the compartment. Each rate for the arrow is always multiplied by the number of individuals in the compartment of origin to calculate the total number of individuals moving out of the compartment at each time step.
For the simple SIR model with population turnover shown above, the equations would be:
It may seem daunting at first to write equations, but it’s actually not so bad once you translate everything into plain English.
| Model term | Plain-language interpretation |
|---|---|
dS/dt |
Change in the number of individuals in the susceptible (S) compartment per unit of time (t).
In this example, time is measured in days.
|
μN |
Number of new individuals entering the susceptible compartment each day, calculated as the birth rate
(μ, individuals per day) multiplied by the total population size (N).
|
-βIS |
Number of individuals leaving the susceptible compartment each day due to infection, equal to the transmission
coefficient (β) multiplied by the number of infectious individuals and the number of susceptible individuals.
|
-αS |
Number of individuals leaving the susceptible compartment each day due to natural death, calculated as the death
rate (α) multiplied by the number of susceptible individuals.
|
You will commonly hear these equations called ordinary differential equations (ODEs) and when people talk about solving ODEs, they largely just mean using the equations to produce graphs showing changes in the number of individuals in each compartment over time.
When you are developing your compartmental model, it is important to have a clear understanding of whether the disease has different impacts or different transmission rates based on the demographic characteristics of the animals. The previous simple SIR model we developed assumed density-dependent transmission dynamics where the likelihood of an animal getting infected increased proportionally to the number of infected individuals. This is probably a reasonable assumption for small populations where your chances of bumping into an infected individual are going to significantly increase the more infected individuals there are.
In large populations, however, this is not usually the case because there is a finite limit to the number of contacts an individual can make. For example, 20,000 people within a single large city of 2 million may have influenza at any given time during winter, but a single individual is unlikely to make 20,000 contacts. Instead, we can keep the contact rate of the individual the same (i.e. 100 total contacts) and just modify the equation so that a proportion of those contacts equal to the disease prevalence (I/N) are actually infectious. For the SIR model with population turnover, the equations would be:
It is important to note the values for β will be different for each system reflecting differences in the assumed contact rate.
The best way to approach this problem is to work compartment by compartment. There should be a + term for every arrow going into a compartment and a – term for every arrow going out of a compartment with the total number of terms for each equation equal to the total number of arrows to/from the compartment. Also remember that if an arrow connects two compartments, what goes out of the first compartment must therefore enter the second compartment.
There are no parameters defined here so let’s make β = transmission coefficient, δ = rate of transition from E to I, α = death rate due to disease, and γ = recovery rate.
The best way to approach this problem is to work equation by equation. There should be a + term for every arrow going into a compartment and a – term for every arrow going out of a compartment with the total number of terms for each equation equal to the total number of arrows to/from the compartment.
In this system, there are two possible outcomes after moving out of the infectious compartment: recovery state or carrier state. The probability of an individual going down either pathway is described by the term θ, which is multiplied by the ϒI term to get the overall rate individuals move into each compartment. Also note that both I and C individuals exert infection pressure on the S compartment.
This is frequency-dependent transmission since we’re not dividing by N
Since the birth rate is equal to the death rate, the population should remain constant over time.
This is describing a situation like Johne’s disease where you can get vertical transmission. To keep the population size constant, we just need to adjust the multipliers for the birth rates going into the S and I compartments.
Read through the description of parvovirus pathogenesis below and develop (1) a compartmental diagram showing the transition of individuals between immunological and infection states and (2) write a series of ODEs to describe the system.
Parvovirus is a highly contagious and common cause of acute, infectious gastrointestinal illness in dogs. Assuming sufficient colostrum ingestion, puppies born to a dam with CPV antibodies are protected from infection for the first few weeks of life; however, susceptibility to infection increases as maternally acquired antibody wanes. Virus is shed in the faeces of infected dogs within 4–5 days of exposure (often before clinical signs develop), throughout the period of illness, and for approximately 10 days after clinical recovery. Infection is acquired through direct oral or nasal contact with virus-containing faeces or indirectly through contact with virus-contaminated fomites (e.g., environment, personnel, equipment). Parvovirus is extremely resistant in the environment and can survive on objects for months to years. Clinical signs of parvoviral enteritis generally develop within 5–7 days of infection but can range from 2–14 days. Initial clinical signs may be nonspecific (e.g., lethargy, anorexia, fever) with progression to vomiting and haemorrhagic small-bowel diarrhoea within 24–48 hr. Physical examination findings can include depression, fever, dehydration, and intestinal loops that are dilated and fluid filled. With appropriate supportive care, 68%–92% of dogs with CPV enteritis will survive. Dogs that recover develop long-term, possibly lifelong immunity.
Parvovirus is a highly contagious and common cause of acute, infectious gastrointestinal illness in dogs. Assuming sufficient colostrum ingestion, puppies born to a dam with CPV antibodies are protected from infection for the first few weeks of life (Animals are born with maternal antibodies: MA compartment); however, susceptibility to infection increases as maternally acquired antibody wanes. (There is a period where animals are susceptible: S compartment) Virus is shed in the faeces of infected dogs within 4–5 days of exposure (often before clinical signs develop) (There is a period where animals have been exposed, but are not yet infectious: E compartment), throughout the period of illness, and for approximately 10 days after clinical recovery (There is a period where animals are actively shedding the virus: I compartment). Infection is acquired through direct oral or nasal contact with virus-containing faeces or indirectly through contact with virus-contaminated fomites (e.g., environment, personnel, equipment). Parvovirus is extremely resistant in the environment and can survive on objects for months to years. Clinical signs of parvoviral enteritis generally develop within 5–7 days of infection but can range from 2–14 days. Initial clinical signs may be nonspecific (e.g., lethargy, anorexia, fever) with progression to vomiting and haemorrhagic small-bowel diarrhoea within 24–48 hr. Physical examination findings can include depression, fever, dehydration, and intestinal loops that are dilated and fluid filled. With appropriate supportive care, 68%–92% of dogs with CPV enteritis will survive(Only a proportion of animals recover and the rest die). Dogs that recover develop long-term, possibly lifelong immunity (Animals recover from the disease and stay immune: R compartment).
What tends to throw people off with this question is having a short period within the infectious period where the animals were showing clinical signs – we only need to model these separately if we are using clinical signs to diagnose disease in the model.
The other things were just making sure that only R animals with antibodies could give birth to MA puppies and accounting for the deaths due to disease.
A deterministic model is one that will produce the exact same result every time you run the simulation because the model parameters, equations, and starting conditions are always the exact same. In biological systems, however, there is a lot of inherent natural variation between individuals in a population that can lead to variability in the epidemic outcomes. For example, the “length of the infectious period” parameter likely has a range of values with some individuals in the population recovering quite rapidly and others taking significantly longer.
A stochastic model is one that incorporates elements of this natural noise either by allowing variability in the time it takes for an event to occur (i.e. sampling from the distribution to determine the length of time an animal spends in a compartment) or allowing different types of events to occur at random (i.e. using Bernoulli trials, which involves assigning all possible outcomes a range of values within 0 to 1 proportional to their likelihood of occurring, generating a random number between 0 and 1, and then choosing the course of action where that value lies). The main reason for including stochasticity is to get an idea for the range of outcomes that could be seen in the real world through the natural variation in biological process. Epidemics often become extinct in small populations due to chance, when the prevalence of infected individuals is low, if there are large variations in the number of susceptible animals, or if the R0 value is low.
In Exercise 2.2, for example, there was a rate ϒ by which individuals left the Infectious compartment and then there was a probability θ that an individual leaving the infectious compartment would become a carrier instead of recovering. If we were programming an individual-based model, we could draw a random number from the distribution of infectious periods for each animal that became infected to determine how long it would spend in the I compartment. Similarly, if the value of θ was 0.275, we would pick a random number between 0 and 1. If the value came back < 0.275, we would move the animal into the Carrier compartment. If the value came back >= 0.275, we would move the animal into the Recovered compartment. These are called Bernoulli trials and you can think about it like flipping a coin to determine the outcome.
It becomes slightly trickier to model stochasticity in traditional compartmental models where the population is aggregated into groups. The Gillespie algorithm is commonly used for this purpose and works by (1) first selecting a random time until the next event occurs, (2) determining what type of event that will be (birth, transmission, recovery, etc), and then (3) moving individuals out of the appropriate compartment at the appropriate proportional rate. It’s easiest to explain this process with a worked example.
Let’s say in our simple SIR model with no population turnover that on day 2 we currently have 100 individuals in the S compartment, 10 individuals in the I compartment, and 20 individuals in the R compartment with β = 0.01 and ϒ = 0.3. The time step for the model is 1 day and there are two events that can occur each day: infection (S -> I) that occurs at the rate βI and recovery (I -> R) that occurs at the rate ϒ. We want to run the model for a total of 50 days to explore the disease dynamics.
From this information, we can predict the number of individuals that will transition between compartments at the next time increment as the result of the events (i.e. numbers moving along each model arrow).
There are 13 total events for the day – 10 individuals that get infected and 3 individuals that recover. The Gillespie algorithm will then:
1. Generate a Poisson process describing the distribution of times until the next event by setting λ = 1/13, plotting out the cumulative distribution function F(x) = 1 – e-λx, and then drawing a time to next event randomly from the distribution. At any given time point, the more individuals that are expected to transition between states, the more rapidly we would expect the next event to happen because the λ value would be smaller. There is a really great tutorial explaining this process here.
In this case, let’s say that the time to next event came back as 0.07 days.
2. Select an event to occur by a similar process to what we saw with the Bernoulli trials. We would expect a 10/13 (76.9%) chance of an S -> I event occurring and a 3/13 (23.1%) chance of an I -> R event occurring. We would normalize this to a 0 to 1 scale and pick a random number between 0 and 1. If the number was < 0.769, the next event would be a S-> I transition and if the number was >= 0.769, the next event would be an I -> R transition. In this case, let’s say the number came back as 0.821 and so we have an I -> R transition.
For disease systems where there are more than two events occurring, you would just divide the range from 0 to 1 into as many compartments as needed and assign an outcome based on where the randomly selected number falls.
3. Update the model by moving the infected individual into the I compartment and increment the time for the model by 0.07 days.
Advancing Animal Care through Shared Learning