Adaptive SIR model with vaccination: simultaneous identification of rates and functions illustrated with COVID-19.

An Adaptive Susceptible-Infected-Removed-Vaccinated (A-SIRV) epidemic model with time-dependent transmission and removal rates is constructed for investigating the dynamics of an epidemic disease such as the COVID-19 pandemic. Real data of COVID-19 spread is used for the simultaneous identification of the unknown time-dependent rates and functions participating in the A-SIRV system. The inverse problem is formulated and solved numerically using the Method of Variational Imbedding, which reduces the inverse problem to a problem for minimizing a properly constructed functional for obtaining the sought values. To illustrate and validate the proposed solution approach, the present study used available public data for several countries with diverse population and vaccination dynamics-the World, Israel, The United States of America, and Japan.

Introduction

Infectious diseases modelling attracted great deal of attention by scientists and people around the world during the COVID-19 pandemic. The novel coronavirus (SARS-CoV-2) started to quickly spread in the early months of 2020 and was announced as a global pandemic by the World Health Organization[1] in March 2020. Population-wide vaccination is critical for achieving herd immunity and for controlling the COVID-19 pandemic while combined with effective testing and preventive measures.

Inevitably, the development of vaccines became the highest priority of governments and pharmaceutical companies[2]. Several vaccines were available in the last months of 2020. As of 30 July 2021, over 28% of the world population is partly or fully vaccinated[3]. Not knowing much about the coronavirus disease in the early months, mathematical models have played an important role in shedding some light on the disease dynamics.

The SIR model categorizes individuals as Susceptible, Infectious, and Recovered. Mathematical infectious disease models based on the classical SIR model[4] are widely used to examine the spread of a disease. These models display compelling results especially during the early period of the pandemic[5-19]. In a recent article[20], the authors studied the immune response of recovered from COVID-19 individuals up to 8 months patients and found that they have considerable immune memory.

Vaccination is a common method of reducing infectious diseases spread[21-25] because it reduces the number of susceptible, from where the reproduction number naturally also decreases. The reproductive number is an indicator of the transmissibility of the virus caused by infectious individuals. It is affected by the population density, vaccinations, quarantines, social distancing, mask wearing and other measures[25,26].

In a recent work[27], the authors present a study of the temporal evolution of epidemic outbreaks accounting for vaccinations with monitored real time COVID-19 data, using the SIRV model, V denoting the relative fractions of currently vaccinated. They make certain assumptions and reduce the time-dependent general SIRV equations to an analytical model. SIR and SEIR models with vaccination are used to simulate and predict the development of the COVID-19 spread, e.g.[25,28-30].

In several very recent publications[31-35] applied to the COVID-19 epidemic, researchers have developed and used SIR and SEIR based models with vaccination to overcome the limitations of the conventional SIR model. The work in[31] presents an investigation of the dynamics of a stochastic SIRV epidemic model with general non-linear incidence and vaccination. The introduced random fluctuations controls the disease outbreak. Zhao et al.[32] use improved SIRV to evaluate the performance of non-pharmaceutical interventions in reducing the number of daily new cases of COVID-19 in South and Southeast Asia. They apply statistical methods to estimate parameters.

The research presented in[33] proposes an SIRV evolutionary game model for infectious disease vaccination strategies based on the scale-free networks with tunable clustering. Their model analyzes the vaccination strategies taking into account factors such as vaccination effectiveness, vaccination cost, treatment cost after illness, government subsidy rate and treatment discount rate. Other notable work is about modelling infectious diseases with herd immunity in a randomly mixed population[34]. The authors formulate two new SIR models to mimic the declining transmission rate of infectious diseases at different stages of transmission. They found that natural herd immunity might not be sufficiently effective in infectious diseases with high reproduction numbers.

Researchers implemented a modification of a SIR model[35] to study the role of the rate of vaccination, rate of transmission and the likelihood of emergence of resistant strains. They used parameters realistically resembling SARS-CoV-2 transmission to run simulations for a total time of three years, with vaccination starting one year into the model.

The goal of the present work is to develop and to demonstrate the effectiveness of an inverse method for identifying the time-dependent functions and parameters of the SIRV model simultaneously (Sect. 2). We apply the method to the adaptive SIRV (A-SIRV) epidemic model using publicly available COVID-19 data. In contrast to other works which use statistical approaches to estimate parameters, we apply an inverse problem approach to identify these parameters. The dynamically estimated rates can be particularly useful in running other simulations, such as in[35], to study the epidemic.

The modified SIR model for the spread of an infectious disease and vaccination (SIRV)

The standard notations in the SIR model are: S(t) denotes the number of susceptible, I(t) – infectives, and R(t) – removed individuals. Assume that the time-dependent function represents the vaccination rate. Then the total number of vaccinated individuals is given by

Every vaccine has different level of efficacy. We assume that the vaccine efficacy impacts disease spread and prevents transmission at the same rate, which is a reasonable assumption according to the study[35]. If people are vaccinated with the same vaccine type / brand, it is possible to introduce the efficacy of the vaccine in the model. Since there is lack of information about the vaccine types and other details concerning the vaccinated individuals, we make the following assumptions, which may not be true for all

Individuals being vaccinated belong to the class S (susceptible) before the vaccine. There are recovered people who vaccinate themselves; this is not included in the model because of lack of data.

Individuals move to the class R (removed) after vaccination. In other words, the model assumes the vaccine is 100% effective against the disease, namely vaccinated people become fully immune.

Due to unavailability of data about COVID-19 variants and vaccine details on a country level, at this stage, the model does not include important assumptions, such as: The transmission rate gives the probability that a random infective person infects a random susceptible person. A major approximation here is the assumption that the population under study is well mixed so that every person has equal probability of coming into contact with every other person. The removal rate gives the probability that an infective person recovers. In the classical SIR model and are constants, while they are functions of time in[14]. The diagram in Fig. 1 describes the adaptive SIR (A-SIR) model with vaccinations (A-SIRV) and time-dependent coefficients:The total population is considered constant in the equations. Fig. 1 shows the diagram for the A-SIRV model, corresponding to the system (2)–(4).

Figure 1

The A-SIRV epidemic model.

Vaccinated individuals could be infected and be infectious if the infection is caused by other variants, known as the vaccine breakthrough problem; there might be multiple variants circulated in the same country, not only a single strain.

Infected individuals could be vaccinated again to improve their immune level; vaccinated individuals could be infected again due to waning of immunity for COVID-19.

The A-SIRV model used in this study assumes that the removed individuals are no longer susceptible nor infectious. The number of cases for recovered from COVID-19 individuals who are re-infected at the present moment is very limited and the rate cannot be estimated; thus, this possibility in not taken into account.

Effective and basic reproductive ratios and

An epidemic occurs if an infective individual introduced into a population of susceptible individuals infects on average more than one other person, namely I(t) is increasing in time. The original SIR model assumes that the transmission and removal rates are constants. Equations (2)–(4), with proper initial conditions, allow the determination of I(t), and S(t), and R(t), if the coefficients and are known constants. However, in the case of a pandemic, the rates may vary in time; hence, and .

The so-called effective reproduction number (also effective reproduction rate or ratio) for a given epidemic is the parameter given byAn epidemic occurs when . Then the fraction of the population that is immune increases (because of vaccination or because of recovering from the disease) so much that , ”herd immunity” is achieved. Hence, the number of new cases occurring in the population will decrease to zero.

Another important characteristic of an epidemic, the basic reproduction number (or ratio, or rate) , is defined as:where N is the size of the total population.

Multiple factors may cause the rates to change over time. In the case of COVID-19, examples include social distancing, restrictions imposed by governments, and preventive treatments. Therefore, we define the effective reproduction number and the basic reproduction number to be functions of time, defined in (5) and (6), respectively.

The main goal of the present work is to identify the time-dependent reproduction rates directly using the SIRV model.

Inverse problem formulation

The initial-value problem consisting of the system of equations (2)–(4), with coefficients and known, along with proper initial conditions derived from the given data, constitutes the direct problem. Note that the vaccination rate function u(t) in equations (2) and (4) can be obtained from the function V(t), which is known from the reported public data on a daily basis.

In reality, the values of the time-dependent parameters and are unknown for a new epidemic disease. Hence, the simultaneous determination of the coefficients and functions from the available data is an inverse problem.

Let [0, P], where P is a number of days, be a time sub-interval (Fig. 2), where approximate values of the susceptible S(t) and currently infectious I(t) are known at specified time moments ,At the same time, the values of are given exactly since the number of vaccinated people is available, as reported in[3]. Therefore, we can approximate the data values of V(t) using a function with a continuous second derivative. From this function, the values can be estimated.

Figure 2

The time nodes where the data values are given.

Since the values of the vaccination rate function u(t) can be found from the available data of vaccinated individuals, we can find I(t), S(t), , and from equations (2), (3), and (7). Note that this problem is inverse and require a special treatment. The approach used here is based on the Method of Variational Imbedding (MVI)[14,36]. Following the idea of MVI, we construct a functional using the original equations (2), (3), and the available data values (7). We define the functional over the sub-interval [0, P] aswhere are the weights prescribed for the i-th node and is the Dirac delta function defined as:

In other words, we substitute the problem for finding the unknown functions I(t), S(t) and the coefficients , in equations (2), (3), and (7) with a problem for minimization of the functional defined by (8). Finding R(t) is straightforward after knowing the values of I(t), S(t), , and .

The absolute minimum of is equal to zero with the functional becoming zero if and only if equations (2), (3) and conditions (7) are satisfied. We want to emphasize here that the functions (S, I, V) and the parameters are unknown and one have to identify them simultaneously by solving a minimization problem.

Method for solving the inverse problem with time-dependent rates

In order to solve the inverse problem with time-dependent rates, we first find the minimum of the functional , defined in (8), numerically over a sub-interval of the entire period assuming constant transmission and recovery rates. The details of solving the minimization sub-problem are given in “Appendix A”.

Let the values of S(t) and I(t) be known at some time moments , shown in Fig. 3, namely

Figure 3

The time nodes and the subsets of fixed length of days for identifying and .

We consider the following two approaches for estimating the coefficients and in the system of equations (2), (3):

Solve the inverse problem for the system of equations (2), (3) under the boundary conditions for . This approach is similar to the method developed in[37] for estimating the coefficient in Euler-Bernoulli equation, later modified in[36] for the SIR equations. It works well if the data represent the exact values of the functions S(t) and I(t). In a real-world data, e.g. the available public data for the COVID-19 pandemic, there usually exists random noise causing oscillations in the numerical results.

Use the solution method, proposed in Sect. 3, for estimating the parameters and as constants on every sub-interval , Then, we use the obtained constant values to approximate the non-constant values and for . This approach is smoothing the data automatically. Knowing the approximate values of the transmission and recovery rates, we obtain the reproduction rates

Results

The numerical simulations have been performed using the available public data from[3] and[38] websites. The number of currently infected persons, , is reported daily on[38]. Both websites[3] and[38], report the total number of infected individuals, from the beginning of the COVID-19 pandemic.

The daily number of vaccinated people, , is taken from[3]. We assume that the reported data is correct. In order to to obtain a smooth approximation of the vaccinated persons V(t) and the vaccination rate , we use cubic spline approximation based on the values . Consequently, the number of susceptible individuals, , for a given country can be found in the following way: , where N is the total population.

We do not pretend that the available data are accurate. For sure, the so-called ”hidden cases” (asymptotic cases or cases without official tests) are not included in the reported data. Moreover, the posted data contain random (human) mistakes which can be treated as random noise. Here, we use the posted data to illustrate the method described in "Method for solving the inverse problem with time-dependent rates" section.

The presented results for the rates , , , and are based on 28-day, 35-day, and 42-day sub-periods, over the entire multi-month period, from August 2020 until August 3, 2021.

The selected countries (Israel, United States, and Japan) for this study represent population with different vaccination dynamics. Israel performed early aggressive vaccination during the first half of 2021, United States is catching up, whereas Japan is behind compared to them. The World is aa good example of aggregated global data. Table 1.

Table 1

Vaccinated, recovered, and susceptible for the World and selected countries as of August 3, 2021.

	Fully vaccinated (%)	Recovered (%)	Susceptible (%)
World	14.78	2.56	82.66
Israel	57.70	9.43	32.87
United States	49.54	10.79	39.67
Japan	30.10	0.75	69.15

The World

Estimated rates , , , and for the World are shown at Fig. 4. According to the reported data as of August 3, 2021: of the population was fully vaccinated; about of the population met the virus. Thus, of the population remained susceptible. The basis and effective reproduction rates have been near 1 since January 2021.

Figure 4

Estimated rates , , , and for the World.

Israel

During the first half of 2021, Israel was the COVID-19 vaccination champion – over of the population of the country was vaccinated by the end of July. At the same time, according to the reports, approximately of the population met the virus; hence, only of the population was susceptible. According to[21], early mass vaccination programs predict a reduction of the effective reproduction rate of infection within communities. The estimated values of the basic and effective reproduction ratios and for Israel are shown in Fig. 5. The effective reproduction rate is practically constant, slightly less then 1 from January to middle June 2021. In June 2021, the rates increased significantly for a short period of time. Then, the reproduction rates decreased but they still remained above 1.

Figure 5

Estimated rates , , , and for Israel.

The United States of America

The reported data for The United States of America state: were fully vaccinated, about of the population met the virus; hence, of the population were susceptible. The estimated values of the transmission rate , , the basic and effective reproduction ratios and are shown at Fig. 6. The transmission rate was the highest in fall 2020. It started to increase again in June-July 2021. COVID-19 is on the rise in many countries, casing a new wave. This surge is due to widespread resumption of normal activities.

Figure 6

Estimated rates , , , and for United States.

Japan

Japan is with about fully vaccinated individuals, approximately reported to be totally infected; this, about are susceptible;. The estimated values of the transmission rate , , the basic and effective reproduction rates and are shown in Fig. 7. While the rates had been relatively reasonable and low until recently, they started to increase lately. This growth can be explained with the 2020 Summer Olympics held from 23 July to 8 August 2021 in Japan.

Figure 7

Estimated rates , , , and for Japan.

Discussion

We performed numerical simulations with the developed method for the A-SIRV model for identifying the transmission, removal and reproduction rates globally and for three countries (Israel, United States, and Japan) using publicly available data. Estimated rates for the World show the presence of three waves and a forth wave being formed since July 2021.

The values for the identified rates for the World do not oscillate, regardless of the number of days P. By contrast, the obtained rates of the selected individual countries have fluctuations during 2021, which are smaller for compared to and . This is because a larger time period is having a smoothing effect of the oscillating functions.

We observe oscillations during time periods that include vaccination data. The vaccination rate u(t) is approximated from the given data for the number of vaccinated individuals V(t). It clearly affects the computed rates in the A-SIRV equations. The oscillations in the rates for Israel, Japan, and to some extent United States can be explained with fluctuations in the vaccination rate u(t) due to irregularities in the data for vaccinated people V(t). We assume that the problem is posed correctly, namely the data have ”physical meaning” and, therefore, a solution of the problem exists, see[39]. The randomness of the vaccination data introduces oscillations in the obtained numerical values of the rates.

Finally, it is important to mention that the A-SIRV model considers only one homogeneous population for the selected countries: thus, giving aggregated results for the estimated rates and the unknown functions. The populations of The World and The United State of America are obviously not homogeneous. Such results can still help see trends of potential future growth of the epidemics and, if needed, guide the design of alternative interventions.

Conclusions

Mathematical models can help with visualizing and predicting the long-term behaviour of an infectious disease, despite of the fact that there are many limitations to using them. For instance, the SIR type models contain many assumptions such as: accuracy and completeness of reported data; mixing of the population; no reinfection; constant population; and so forth. The present work studies the performance of a method for an epidemic based on an inverse problem approach for estimating the time-dependent transmission and removal rates in the A-SIRV epidemic model. The inverse problem is solved by defining a minimization problem using the entire dataset for the examined population, with available COVID-19 data. The work utilizes an inverse problem approach to the time-dependent transmission and removal rates identification as well as the unknown functions in the A-SIRV (Adaptive SIRV) model. This can give insight into how well the method identifies the parameters that can be used to predict the infectious disease spread. If conditions change, then the predictions may no longer be accurate; hence, adjustments will be required based on the existing conditions, for obtaining new predictions.

Table 2

Obtained values of the constants and and the rate of convergence for four different values of the mesh spacing .

τ	Nβ	rate_β	γ	rate_γ
1/10	0.1005883672528720	–	0.0741311130168358	–
1/20	0.1005886480876050	–	0.0741311570649426	–
1/40	0.1005887182822550	2.0002884	0.0741311680719028	2.0006639
1/80	0.1005887358165650	2.0011804	0.0741311708184778	2.0027105
1/160	0.1005887402874910	1.9715351	0.0741311715365409	1.9354511