SIR Model

  • S: Susceptible is a group of people who are vulnerable to the infection.
  • I: Infectious population is the group of people who can transmit the infection.
  • R: Recovered(and dead) population is the group of people who are no longer infectious, this includes both the recovered population and people succumbed to the epidemic.

A key assumption of the model is a recovered person generates immunity towards the epidemic.

The model also assumes that the population remains constant, i.e, there are no new births and no deaths on account of reasons other than the epidemic.

Initial Conditions

Total population(N) is assumed to be 1.

Initial state for infectious category(I0) is the fraction of total population which is infectious at time T0.

Initial state for susceptible population(S0) is the remaining population(N-I0) assuming that no person has been vaccinated.

It is also assumed that there is no recovered person in the beginning.

Epidemic starts to spread

An infection turns into an epidemic if R0 > 1, it dies out if R0 < 1.

In [1]:
#import libraries
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

import warnings
warnings.filterwarnings("ignore")
In [2]:
#Function to model the SIR framework
def sir_model(I0=0.01, beta=0.6, gamma=0.1):
    """
    Function will take in initial state for infected population,
    Transmission rate (beta) and recovery rate(gamma) as input.
    
    """
    N=1          #Total population
    I=I0         #Initial state of I default value 1% of population
    S=N-I        #Initial state of S
    R=0          #Initial State of R
    C=I          #Initial State of Total Cases
    beta=beta    #Transmission Rate
    gamma=gamma  #Recovery Rate

    inf=[]       # List of Infectious population for each day
    day=[]       # Time period in day
    suc=[]       # List of Susceptible population for each day
    rec=[]       # List of Recovered population for each day
    conf=[]      # List of Total Cases population for each day
    
    for i in range(60):
        day.append(i)
        inf.append(I)
        suc.append(S)
        rec.append(R)
        conf.append(C)

        new_inf= I*S*beta       #New infections equation (1)   
        new_rec= I*gamma        #New Recoveries equation (2)
        
        I=I+new_inf-new_rec     #Total infectious population for next day
        S=S-new_inf             #Total infectious population for next day
        R=R+new_rec             #Total recovered population for next day
        C=C+new_inf             #Total confirmed cases for next day
    
    max_inf=round(np.array(inf).max()*100,2)     #Peak infectious population in percentage
    max_conf=round(np.array(conf).max()*100,2)   #Overal infected population in percentage
    
    print(f"Maximum Infectious population at a time :{max_inf}%")
    print(f"Total Infected population :{max_conf}%")
    
    #Visualizing the model
    sns.set(style="darkgrid")
    plt.figure(figsize=(10,6))
    plt.title(f"SIR Model: R0 = {round(beta/gamma,2)}", fontsize=18)
    sns.lineplot(day,inf, label="Infected")
    sns.lineplot(day,suc,label="Succeptible")
    sns.lineplot(day,rec, label="Recovered")
    #sns.lineplot(day,conf, label="Confirmed") #Generally total infected population is not plotted 
    plt.legend()
    plt.xlabel("Time (in days)")
    plt.ylabel("Fraction of Population")
    plt.show()
In [3]:
#Plot infectious population and total infected population for various R0 values
def sir_model_betalist(I0=0.01, betalist=[0.5,0.8], gammalist=[0.15,0.25,0.5]):
    """
    Function takes Initial Infected Population(I0), list of transmission rates (betalist)
    and list of recovery rates(gammalist) as arguments.
    Plots Infectious population and Infected Population vs time for input parameters
    """
    
    for gamma in gammalist:
        # Plotting Infectious Population
        plt.figure(figsize=(10,6))
        sns.set(style="darkgrid")
        plt.title("SIR Model: Infectious Population", fontsize=18)
        
        for beta in betalist:
            N=1
            I=I0
            S=N-I
            gamma=gamma
            R0=beta/gamma
            
            inf=[]
            day=[]
            for i in range(50):
                day.append(i)
                inf.append(I)
                new_inf= I*S*beta
                new_rec= I*gamma
                I=I+new_inf-new_rec
                S=S-new_inf
            
            inf_max=round(np.array(inf).max()*100,1)
            sns.lineplot(day,inf, label=f"R0: {round(R0,2)} Peak: {inf_max}%")
            plt.legend()
        plt.show()
        
        # Plotting Total Infected Population
        plt.figure(figsize=(10,6))
        plt.title("SIR Model: Total Confirmed Cases", fontsize=18)       
        for beta in betalist:
            N=1
            I=I0
            S=N-I
            C=I
            gamma=gamma
            R0=beta/gamma
            day=[]
            conf=[]
            for i in range(50):
                day.append(i)
                conf.append(C)

                new_inf= I*S*beta
                new_rec= I*gamma
                I=I+new_inf-new_rec
                S=S-new_inf
                C=C+new_inf
            conf_max=round(np.array(conf).max()*100,1)
            sns.lineplot(day,conf, label=f"R0: {round(R0,2)} Total :{conf_max}%")
            plt.legend()
        plt.show()

According to this, R0 is between 2.2 and 5.7. So I will model a scenario where R0 is 6, followed by 3, and with an initial state of infectious population at 1%.

In [4]:
sir_model(beta=0.9,gamma=0.15)
Maximum Infectious population at a time :59.6%
Total Infected population :99.94%
In [5]:
sir_model(beta=0.45, gamma=0.15)
Maximum Infectious population at a time :32.22%
Total Infected population :95.01%
In [6]:
sir_model(beta=0.3,gamma=0.15)
Maximum Infectious population at a time :16.42%
Total Infected population :79.05%
In [11]:
sir_model(beta=0.45,gamma=0.3)
Maximum Infectious population at a time :7.28%
Total Infected population :60.07%

For NY, according to rt.live, the R0 is 0.83.

"Flattening the Curve"

Flatten the Curve is a phrase which is very familier in these day, this is what we have done above. It can be observed that, by decreasing R0, we reduce the spread of the infection and the peak level of infected population. Decreasing R0 will also delay the time for peak level of infection, buying us valuable time to prepare tackle the infection effectively.

Let's now check how the peak level of infection and the extent of infection varies with R0.

In [10]:
sir_model_betalist(I0=0.01,betalist=[0.166,0.2, 0.25,0.30,0.4,0.6,0.8,1,1.2], gammalist=[0.20])

As R0 decreases, the curve for infectious population flattens out. When R0 is 6, more than 60% of the population is infections by around 8th day, whereas when R0 is 2, maximum infectious population is 16.7% by around 24th day. That mean reduction in R0 from 6 to 2 gives us additional two weeks to prepare to tackle the peak infection. It can also be observed that the extent of infection also decreases with reduction in R0. When R0=6, total population gets infected whereas when R0=2, 20% of the population remains uninfected when the epidemic is over.

Despite COVID-19 being a relatively less fatal infection, it wreaked havoc all over the world as it quickly spread and caught many countries unprepared, overwhelming their health infrastructure. Health infrastructure is said to be overwhelmed when the number of infectious population is above the total beds available at the hospitals. Yes, the health infrastructure could be ramped up rapidly by converting other infrastructure to makeshift hospitals and calling retired health professionals back into duty etc. as followed by many countries but it won’t be able to keep up with if the cases increase exponentially. The faster the infection curve rises, the quicker the local health care system gets overloaded beyond its capacity to treat people. As we're seeing in Italy, Spain, United States etc. more and more new patients may be forced to go without ICU beds, and more and more hospitals may run out of the basic supplies they need to respond to the outbreak.

How to Flatten the Curve

The curve can be flattened out by reducing R0 which can be achieved by:

  • Reducing susceptible population (S0),
  • Reducing transmission rate (beta)
  • Increasing recovery rate (gamma)


Reducing Susceptible population
The best method for decreasing susceptible population is vaccination, however no vaccination is available for COVID-19. Overall susceptible population can be contained in a region level by implementing travel restrictions, lock-downs etc.
Increasing recovery rate (gamma)

Increasing recovery rate depends on health infrastructure and treatment methods. Significant improvement in this area is both time consuming and capital intensive.
Reduction in transmission rate

Transmission rate depends on the number of persons from susceptible population an infected individual comes into contact (contact rate, k) and probability that a contact becomes infected (transmissibility, t).
Reduce the contact rate by social distancing measures like closing of schools, offices etc., suspension of public transport, self-isolation of susceptible individuals, quarantine etc.
Reduce the transmissibility by improving personal and environmental hygiene by frequent hand washing, covering nose and mouth while sneezing, disinfecting environment, use of face mask etc.

In [ ]: