Time-continuous and time-discrete SIR models revisited: theory and applications.

Since Kermack and McKendrick have introduced their famous epidemiological SIR model in 1927, mathematical epidemiology has grown as an interdisciplinary research discipline including knowledge from biology, computer science, or mathematics. Due to current threatening epidemics such as COVID-19, this interest is continuously rising. As our main goal, we establish an implicit time-discrete SIR (susceptible people-infectious people-recovered people) model. For this purpose, we first introduce its continuous variant with time-varying transmission and recovery rates and, as our first contribution, discuss thoroughly its properties. With respect to these results, we develop different possible time-discrete SIR models, we derive our implicit time-discrete SIR model in contrast to many other works which mainly investigate explicit time-discrete schemes and, as our main contribution, show unique solvability and further desirable properties compared to its continuous version. We thoroughly show that many of the desired properties of the time-continuous case are still valid in the time-discrete implicit case. Especially, we prove an upper error bound for our time-discrete implicit numerical scheme. Finally, we apply our proposed time-discrete SIR model to currently available data regarding the spread of COVID-19 in Germany and Iran.

Introduction

Motivation

Since its outbreak in Wuhan (China) in December 2019, the quick spread of COVID-19 has threatened people worldwide. Politicians around the globe have to balance different interests and need to make tremendous decisions which impact our daily life. For these reasons, governments around the world heavily rely on scientific advice in the present situation. Thus, John Hopkins University has been collecting epidemiological data regarding COVID-19 from many countries during the last months [1, 2]. Additionally, many biological and medical studies regarding different aspects of this new coronavirus have been rapidly appearing in scientific journals [3-9]. However, to estimate the impact of COVID-19, governments need forecasts from as many models as possible.

Kermack and McKendrick introduced the now well-known SIR model in one of mathematical epidemiology’s most groundbreaking works in 1927 [10]. They assumed a fixed population size and divided this population into three different homogeneous groups of people, namely susceptible people, infectious people, and recovered people, excluding births, deaths, and deaths by disease from their model. Due to its success and simplicity, their work was reprinted in 1991 [11-13]. In upcoming decades, epidemiologists and mathematicians have developed many variants and extensions of this basic model by, for example, adding age or spatial structures [14-20].

After the outbreak of COVID-19, many scientists have recently published articles with emphasis on epidemic forecasts which strongly relate to mathematical models. Many approaches, mainly focusing on stochastic arguments, with respect to predicting forecasts of the total number of infected people have appeared during the last weeks [21-26] or in the past [27, 28]. Recently, neural networks have been applied to forecasting [29] and, additionally, different authors such as Atangana, Baleanu, or Khan have used fractional differential equations in extended SIR-type models to investigate the spread of COVID-19 or mathematical biology in general [30-34].

Contributions

There are works with respect to SIR models [35-37] and their time-discrete versions [38]. However, one finds mainly explicit schemes with respect to time-discrete SIR models in the aforementioned works and references therein. For this reason, we summarize and extend some results on properties of the time-continuous classical SIR model, and we propose an implicit time-discrete version of this classical SIR model and prove that this time-discrete variant maintains many of time-continuous version’s properties. Hence, the aim of our study is to propose a nonautonomous SIR model, analyze the properties of its time-continuous formulation, and develop an implicit numerical solution algorithm where the main properties of the time-continuous variant are conserved. More precisely, our main contributions can be summarized as follows:

At first, we propose a modified time-continuous SIR model with time-varying transmission and recovery rates;

Secondly, we conclude well-posedness of our time-continuous problem formulation. This includes global existence in time, global uniqueness in time, based on an inductive application of Banach’s fixed point theorem, and continuous dependence on initial conditions and time-varying rates;

In case of the time-discrete implicit model, we provide unique solvability, monotonicity properties, and an upper error bound between the solution of the implicit time-discrete problem formulation and the solution of the time-continuous problem formulation;

Finally, we compare our model predictions with real-world data regarding the spread of COVID-19 from different countries. Data have been collected by John Hopkins University.

Conclusively, we summarize our results and give a brief outlook on possible extensions.

Time-continuous SIR model

In this section, we portray the time-continuous SIR model and its properties.

Mathematical background material

Here, we recall Lipschitz continuity of a function on Euclidean spaces.

Definition 1

([39, Sect. 3.2])

Let . If , a defined function is called Lipschitz continuous on S if there exists a nonnegative constant such that holds for all . Here, denotes a suitable norm on the corresponding Euclidean space.

Let be open, let . We shall call F locally Lipschitz continuous if for every point there exists a neighborhood V of such that the restriction of F to V is Lipschitz continuous on V.

In a more general framework, we consider a nonlinear initial value problem where we define our solution vector , our vectorial function , and our initial condition . To conclude global existence in time, we can apply the following theorem that is a direct consequence of Grönwall’s lemma.

Theorem 1

([39, Theorem 4.2.1])

If is locally Lipschitz continuous and if there exist nonnegative real constants B and K such that holds for all , then the solution of the initial value problem (2) exists for all time and, moreover, it holds for all .

To prove global uniqueness in time, we need Banach’s fixed point theorem since fixed point theorems have been successfully applied as a versatile tool in different areas of mathematics [40, 41].

Theorem 2

([42, Theorem V.18])

Let be a complete metric space with the metric mapping . Let be a strict contraction, i.e., there exists a constant such that holds for all . Then the mapping T has a unique fixed point.

Finally, since we want to provide continuous dependence on initial conditions and other data, we need the following inequality named after Gronwall and Bellman.

Theorem 3

([43, Theorem 1.3.1])

Let . Let u, be continuous and nonnegative functions. Let be a continuous, positive, and nondecreasing function. If the inequality holds for all , then we obtain for all .

Continuous problem formulation

At first, let us state the model’s assumptions [16, 17]: We briefly comment on our choice of time-varying transmission rates. The choice of time-dependent transmission rates is possible because countermeasures such as lock-downs, social distancing, or other political actions such as curfews reduce possible contacts between susceptible and infectious people. In addition to that, time-varying recovery rates might also be an interesting choice due to different medical treatments over the course of new epidemics such as COVID-19.

The population’s size N is fixed over time, i.e., for all ;

We divide a population into three homogeneous subgroups, namely susceptible people (S), infectious people (I), and recovered people (R). We can clearly assign every individual to exactly one subgroup. Hence, we obtain for all ;

Additionally, no births and deaths occur;

The time-varying transmission rate is Lipschitz continuous and continuously differentiable once, and there exist positive constants and such that holds for all ;

The time-varying recovery rate is Lipschitz continuous and continuously differentiable once, and there are positive constants and such that holds for all .

Furthermore, we exclude age structures or spatial relationships from our time-continuous model [16, 19]. For abbreviation, we write . Our equations of the time-continuous SIR model read as follows: with initial conditions , , and . We portray a chart of the flow between the different three subgroups in Fig. 1. Obviously, the equation is valid, and hence the first assumption is fulfilled.

Figure 1

Flowchart of the three different subgroups described by the time-continuous SIR model

Nonnegativity and boundedness of solution

Now, we prove boundedness of the solution. For this purpose, we modify ideas from [17] and [44] because we, in contrast to them, consider bounded, time-varying transmission rates and recovery rates .

Lemma 1

Each solution function of (8) is bounded below by zero.

Proof

Here, we split the proof into three parts.

1) We first consider . Separation of variables leads to Integration yields and this implies Hence, it holds for all .

2) Let us continue with . Separation of variables gives us and, applying the same argument as in our first step, we conclude This yields for all .

3) Finally, since holds, we clearly obtain and for all is valid.

This completes our proof. □

Since holds by our first assumption, we can finally state our boundedness theorem.

Theorem 4

For all solution functions of (8), we have: for all .

;

;

Global existence in time

In contrast to many other works, we formulate a theorem regarding global existence in time of (8) based on Theorem 1. For abbreviation, we introduce the supremum norm on an arbitrary time interval for a continuous function . A similar definition holds for vector-valued functions. In our case, we consider in general. Now, we show that a possible solution to (8) exists for all times .

Theorem 5

The nonlinear first order ordinary differential equation system (8) has at least one solution which exists for all .

By defining , we can set Clearly, G is Lipschitz continuous. By considering the supremum norm on our Euclidean space and applying the triangle inequality, we get from (9) by the boundedness of our solution functions and the boundedness of our time-varying transmission and recovery rates. Thus, all our assumptions of Theorem 1 are fulfilled and our proof is complete. □

Global uniqueness in time

We present our proof of global uniqueness in time based on an inductive application of Banach’s fixed point theorem, i.e., that the initial value problem (8) is uniquely solvable for all times .

Theorem 6

The nonlinear first order ordinary differential equation system (8) has a unique solution that exists for all .

1) Let us first consider the time interval where we have to choose τ accordingly such that Banach’s fixed point theorem is applicable.

2) We need one brief lemma. Let be arbitrary. By zero addition and application of the triangle inequality, we obtain

3) We assume that S, I, R, S̃, Ĩ, are two solutions of (8). At first, it holds by our inequality of the second step.

4) Secondly, we obtain by our inequality of the second step and application of the triangle inequality.

5) Furthermore, we conclude that holds.

6) Summarizing the previous steps, we obtain If we choose , this implies and hence we conclude the uniqueness of solution on the time interval .

7) Inductively, we see that we can derive this contraction property on all time intervals for all by choosing as our starting point and for our initial conditions. Henceforth, our proof of global uniqueness in time is complete. □

Continuous dependence on initial conditions and time-varying rates

Here, we consider the perturbed initial value problems with initial conditions , , and with initial conditions , , , where , , , are different time-varying transmission and recovery rates. Now, we prove that small perturbations in initial conditions or small differences in time-varying transmission or recovery rates lead to small differences in the solutions on short time intervals . This fact is summarized in the following theorem.

Theorem 7

Let and be the solutions of (10) and (11). Define the function and the constant We see that holds for arbitrary with given .

1) Let us first mention that we often use the inequality for arbitrary as proven in Theorem 6. Additionally, we see that hold for all .

2) At first, we obtain the inequality for arbitrary by application of the triangle inequality, boundedness of our all functions, and the inequality of our first step.

3) Secondly, we have to estimate . We see that holds for arbitrary . The summand I can be estimated in the previous step. This yields for arbitrary . For the third summand II, we observe that is valid for arbitrary t. Summarizing these results, we obtain for arbitrary .

4) Now, we must estimate . It holds for arbitrary .

5) Finally, we obtain for arbitrary . This implies

6) Define the functions Since all the assumptions of Theorem 3 are fulfilled, we see that holds for arbitrary , which finishes our proof. □

Monotonicity properties and long-time behavior

We now investigate the long-time behavior of solution, some monotonicity properties and summarize our results in the following theorem.

Theorem 8

We get: for all solution functions of (8).

S is monotonically decreasing, and there exists a number such that . It even holds ;

R is monotonically increasing, and there exists a number such that ;

I is Lebesgue-integrable on and

1) Since for all and by Theorem 4, is monotonically decreasing and bounded below by zero. This implies the existence of such that . Additionally, by considering for , we obtain and separation of variables implies Integration yields for all . Hence, it holds .

2) Since for all and is true by application of Theorem 4, is monotonically increasing and bounded above by N. This yields the existence of such that holds.

3) Since holds, integration yields because all functions α, S, are bounded and nonnegative. Therefore, we obtain that I is Lebesgue-integrable on and . This finishes our proof. □

Calculation of the time-continuous basic reproduction number

In our nonautonomous SIR model, the time-dependent basic reproduction number can be defined by which is similar to the constant case [17, 45, 46].

Lemma 2

Equation (13) is well defined.

We observe that is valid for all . This proves our claim. □

Time-discrete implicit SIR model

In this section, we examine time-discrete versions of the given time-continuous SIR model (8). Let us assume that our time interval can be divided by a strictly increasing sequence for with and . For abbreviation, we write for all and an arbitrary time-dependent function f.

Discussion of formulations

Here, we only state a fully explicit scheme and a fully implicit scheme of the time-continuous SIR model (8) for all . Both formulations fulfill for all . However, the fully explicit scheme (14) simply reduces to a linear system, while the fully implicit scheme (15) maintains the nonlinear structure of the continuous problem formulation (8). For this reason, we investigate this fully implicit scheme in the following.

Time-discrete implicit problem formulation

We assume that and are given for all and that for all and that , , and are given. As later observed in our numerical examples, these assumptions are fulfilled in epidemiological data of the spread of COVID-19. An implicit solution scheme of (15) reads as follows: for all . Now, we are able to obtain an appropriate solution scheme from (17) which even implies unique solvability for all under the assumption that , , and for all .

Unique solvability

Our main ingredient is the equation from (17). Plugging into (18) and writing yields for all . Hence, we get and by setting and we get and can finally conclude for all . We now have an explicit solution formula for for all and therefore also for and for all . Summarizing our results, we can formulate the following theorem.

Theorem 9

Let us assume that and are given for all , that holds for all , and that , , and are prescribed. The implicit solution scheme (17) is uniquely solvable for all . It holds (22) for all with and from (20) and (21).

We show that many of the continuous properties from Theorems 4 and 8 even translate to the time-discrete implicit scheme (17).

Theorem 10

For our time-discrete implicit solution scheme (17), we have:

for all ;

for all and for all ;

for all and for all ;

.

1) It holds due to (22) and due to (16) for all .

2) By our first property and due to (16), we have the inequality for all . Again by our first property, we obtain for all .

3) By our first property and due to (16), we obtain the inequality for all . Again by our first property, we conclude for all .

4) Since is monotonically increasing and bounded above by the total population size N, there exists a nonnegative constant such that . Furthermore, it holds which yields This implies and completes our assertion’s proof. □

Error analysis

Now, we want to provide an upper bound for error propagation. Before proving this statement, we need to formulate some assumptions for our convergence analysis. We summarize them in the following list:

Let be the considered time interval where ;

Let the initial conditions of the time-continuous and the time-discrete models coincide;

Let the solution functions S, I, be twice continuously differentiable;

Let the time-varying transmission rate and the time-varying recovery rate be once continuously differentiable;

Let the time-varying transmission and recovery rates be bounded, i.e., there are nonnegative constants , , , such that and hold for all ;

Choose for all and set .

Under these conditions, we obtain the following theorem where we adapt ideas from the error analysis of an explicit-implicit solution algorithm as presented in [20].

Theorem 11

If the aforementioned assumptions are fulfilled, the difference between the solution of the time-continuous problem formulation (8) and the solution of the time-discrete problem formulation (17) fulfills

We briefly describe our strategy first because this proof is technical. We begin with an estimation of local errors between time-continuous and time-discrete solutions. Afterwards, we consider error propagation in time. Conclusively, we investigate the cumulation of these errors. Time-discrete solutions are written as at time and time-continuous solutions as at the same time.

1) For examination of local errors, we assume that hold for arbitrary on the time interval . Here, we consider solely one time step and denote corresponding time-discrete solutions by , , and .

1.1) It first holds and this implies by application of the triangle inequality. We estimate these terms separately. For , we obtain and the mean value theorem of calculus implies the existence of such that holds. This yields For , we see that is valid by the definition of , boundedness of our solution functions, and application of the triangle inequality. By plugging into , we obtain by the boundedness of our solution functions and application of the triangle inequality. By the mean value theorem, there exists such that holds. Hence, an additional application of the triangle inequality and boundedness of the solution functions yields Plugging (25) into , we conclude that holds. Combining (24) and (26), we infer that is valid. Thus, it holds

1.2) We observe that is valid. Hence, it follows by application of the triangle inequality and with similar arguments as provided in the previous step. As holds, we further obtain by the boundedness of the solution functions and application of the triangle inequality. By using we obtain by the boundedness of the solution functions, the mean value theorem of calculus, and application of the triangle inequality. Additionally, it holds by application of the mean value theorem of calculus. Combining (28) and (29) and plugging these results into yields Plugging this inequality into implies

1.3) We see that holds. This implies By applying and we obtain Plugging this inequality into yields

1.4) Define . It holds for local errors on time intervals .

2) In reality, , , and do not exactly lie on the graph of the time-continuous solution. Therefore, we must examine how procedural errors such as , or propagate to the th time step. These investigations are carried out in step 2) and in step 3). By (9), we see that holds, and this implies We see that holds, and this implies Hence, we conclude with .

3) Now, we want to prove the upper error bound between the time-discrete solution and the time-continuous solution. At first, we notice that is valid for by (32), by (33), and by our assumption that initial conditions of the time-continuous and the time-discrete models coincide. For , we obtain For arbitrary , we assume that is valid. This yields by induction. By applying the geometric series, we obtain If we assume , we conclude and hence it follows which finishes our proof of (23). □

Calculation of the time-discrete basic reproduction number

In our nonautonomous time-discrete SIR model, the time-dependent basic reproduction number can be defined by for arbitrary , which is similar to the case of constant transmission and recovery rates [17, 45].

Lemma 3

Equation (35) is well defined.

This proof is identical to Lemma 2. □

Numerical algorithm

We are now able to give a brief description of our numerical algorithm to solve the time-discrete implicit solution scheme (17). Here, we summarize our inputs, our computational steps, and our algorithmic outputs. We sketch the resulting algorithm in Table 1.

Table 1

Numerical algorithm for the time-discrete implicit SIR solution scheme (17)

Inputs:	– Population size N
	– Initial values S_1>0, I_1>0, and R_1≥0
	– Time-varying transmission rate sequence {α_j}j=2M
	– Time-varying recovery rate sequence {β_j}j=2M
	– Strictly increasing sequence {t_j}j=1M of time points with t_1=0 and t_M=T
Step 1:	– Compute all Δ_j+1=t_j+1−t_j for all j∈{1,…,M − 1}
Step 2:	– Compute I_j+1 by (22), (21) and (20) for all j∈{1,…,M − 1}
	– Compute S_j+1 and R_j+1 by (17) for all j∈{1,…,M − 1}
Outputs:	– Sequences {S_j}j=1M, {I_j}j=1M and {R_j}j=1M

Numerical examples with discussion

We apply our time-discrete implicit SIR solution scheme (22) from Table 1 to available data regarding the spread of COVID-19 in Germany and Iran from John Hopkins University [1, 2]. These countries are chosen because they update confirmed, dead, and estimated recovered cases on a regular basis. In Table 2, we summarize projected population sizes for 2019 from the United Nations [47].

Table 2

Projected population sizes for Germany and Iran for 2019

Country	Germany	Iran
Population size	83,784,000	83,993,000

At first, we consider the example of Germany in detail. We thoroughly describe our approach to check our model’s validity. We only give some simple parameter estimation techniques and vary user-chosen time-dependent parameter functions. Since inverse problems are an active field of research, we refer the readers to works by Bock and Schittkowski for more sophisticated parameter estimation techniques in dynamical systems [48, 49]. Our work also focuses on usefulness in possibly describing real-world data. Afterwards, we state some computational results for data from Iran.

Description of our approach by the example of Germany

Data preprocessing

To apply our model, we have to process the given data of cumulative confirmed infected people , cumulative confirmed dead people , and cumulative confirmed recovered people . For our model, we need to compute the processed real-world data for all . The unprocessed and processed data for Germany are depicted in Figs. 2 and 3. On the one hand, it can be clearly seen in Fig. 2 that both sequences of cumulative infected and cumulative recovered people are monotonically increasing for our unprocessed data. On the other hand, we notice in Fig. 3 that the behavior changes for our processed data.

Figure 2

Unprocessed data for Germany with (1 March 2020) and (1 September 2020)

Figure 3

Processed data for Germany with (1 March 2020) and (1 September 2020)

Calculation of time-varying transmission and recovery rates from real-world data

Here, we present an algorithm to calculate our time-varying transmission and recovery rates based on our numerical algorithm (15) for all . We rely on the first equation and the last equation for susceptible and recovered people. Short calculations with the assumptions and yield and Summarizing our results, we obtain and and a short algorithmic summary can be found in Table 3.

Table 3

Numerical algorithm for the time-discrete implicit SIR solution scheme (17)

Inputs:	– Population size N
	– Real-world data {S_j˜˜}j=1M, {I_j˜˜}j=1M, and {R_j˜˜}j=1M according to (36)
	– Strictly increasing sequence {t_j}j=1M of time points with t_1=0 and t_M=T
Step 1:	– Compute all Δ_j+1=t_j+1−t_j for all j∈{1,…,M − 1}
Step 2:	– For all j∈{2,…,M}, compute α_j˜˜ and β_j˜˜ according to (37) and (38) with real-world data
Outputs:	– Sequences {α_j˜˜}j=2M and {β_j˜˜}j=2M

The time-varying transmission rate from real-world data for Germany is presented in Fig. 4. Clearly, the transmission rate decreases due to countermeasures such as local lock-downs and voluntary social distancing by the population. However, a weak increasing trend can be seen at the end of the time-series in August. Possible explanations might be opening of schools, universities or people who do not wear masks for protection.

Figure 4

Time-varying transmission rate from real-world data for Germany with (1 March 2020) and (1 September 2020)

The time-varying recovery rate from real-world data for Germany is depicted in Fig. 5. At the early stage of an epidemic, there are possible just few recoveries, thus this rate is relatively small. After some time, this situation changes as more people defeat the disease and recover. The rate seems to be constant with heavy variations due to the test capacity. Additionally, there are unknown cases because these people might have a mild disease course.

Figure 5

Time-varying recovery rate from real-world data for Germany with (1 March 2020) and (1 September 2020)

Calculation of time-dependent basic reproduction number from real-world data

Now, the time-dependent basic reproduction number is readily computed by (35). Our computational results from this approach are portrayed in Fig. 6. Since there are only few recovered people at the beginning of disease, our computations provide high numerical basic reproduction number. We observe that the computational basic reproduction number is monotonically decreasing in spring due to political countermeasures and social distancing. However, the graph shows that the computation basic reproduction number rises in summer because contacts between people rose. Variations are seen for similar reasons as mentioned for our time-dependent transmission and recovery rates.

Figure 6

Time-varying basic reproduction number from real-world data for Germany with (1 March 2020) and (1 September 2020)

Parameter estimation—a simple least-squares approach

We only briefly sketch the parameter identification problem because it is an inverse problem [50, 51]. A deep discussion is beyond the scope of this paper and would be a topic of own interest. By looking at Figs. 4 and 5, we assume and with real constants , , and β which we determine from real-world transmission and recovery rate sequences and . Since and for all , we can assume and . Since is nonlinear, we use the transformation with and as in the case of maximum log-likelihood estimation. Now, our cost function reads as follows: We obtain the following theorem.

Theorem 12

Assume that holds for the strictly increasing time sequence . The cost function (42) possesses a unique local minimizer. In fact, this unique local minimizer is even a unique global minimizer.

Proof

1) We first show that J possesses a unique local minimizer.

1.1) To achieve our goal, we calculate the first derivatives. We obtain

1.2) To get a local minimizer , it is necessary that all partial derivatives vanish at candidates for local extrema.

1.2a) Setting , we conclude

1.2b) Setting , we infer that holds.

1.2c) If we set , we obtain and this yields Finally, we conclude that holds.

1.3) Since the Hessian is given by we investigate all determinants of all upper-left sub-matrices. The first determinant because holds. The second determinant reads as follows: By the Hölder inequality, we obtain with equality only in the case that all are equal. Since we have a strictly increasing time sequence, the second determinant of the upper-left sub-matrices is also positive. Finally, the determinant of the full matrix is positive as well. Hence, our cost function J is strictly convex. Conclusively, it possesses a unique local minimizer by [52, Theorem 2.4].

2) By strict convexity, we infer that the unique local minimizer is also the unique global minimizer of our cost function J by [52, Theorem 2.5]. □

We summarize our algorithmic approach for parameter estimation of our time-varying transmission and recovery rates in Table 4.

Table 4

Numerical algorithm for our parameter estimation approach of our cost function (42)

Inputs:	– Strictly increasing sequence {t_j}j=1M of time points with t_1=0 and t_M=T
	– Sequences {α_j˜˜}j=2M and {β_j˜˜}j=2M
Step 1:	– Compute β̂ by (43)
Step 2:	– Compute γ_2ˆ by (45)
Step 3:	– Compute γ_1ˆ by (44)
Step 4:	– Compute α_1ˆ and α_2ˆ according to transformation (41)
Outputs:	– Parameters α_1ˆ, α_2ˆ, and β̂ for our parametric rates (39) and (40)

Results for our parameter estimation approach using German data for short-term predictions

In Fig. 7, we see that the assumption of exponentially decaying time-dependent transmission rates is acceptable at the beginning of spreading disease with respect to German data. Due to short-term prediction, we notice that the constant recovery rate is underestimated in Fig. 8. Conclusively, both assumptions seem to be acceptable at the first weeks of a spreading disease. Computational results for two models on the time interval are depicted in Figs. 9–12. Figures 9–12 indicate that sensitivity of parameters is really an issue in epidemiological models. This is in accordance with Theorem 7. These results also imply that an exponentially decaying transmission rate is an acceptable choice at the beginning of spreading disease.

Figure 7

Time-varying transmission rates from real-world data and from parameter estimation for short-term prediction of German data with (1 March 2020) and (2 May 2020). The estimated parameters for model 1 are and . For model 2, we fix and choose

Figure 8

Recovery rates from real-world data and from parameter estimation for short-term prediction of German data with (1 March 2020) and (2 May 2020). The first estimated recovery rate reads for the mean value on the full interval. The second estimated recovery rate reads as the mean value on the time interval because the fluctuations in β arise while there is no recovery rate estimation possible for the first days

Figure 9

Computational results for model 1 of infected people in Germany for , , for short-time simulation on (26 March 2020) and (2 May 2020) with real-world data as initial conditions

Figure 12

Computational results for model 2 of recovered people in Germany for , , for short-time simulation on (26 March 2020) and (2 May 2020) with real-world data as initial conditions

Computational results for model 1 of recovered people in Germany for , , for short-time simulation on (26 March 2020) and (2 May 2020) with real-world data as initial conditions

Computational results for model 2 of infected people in Germany for , , for short-time simulation on (26 March 2020) and (2 May 2020) with real-world data as initial conditions

Results for our parametric approach using German data from May to September

Figures 4 and 5 indicate that constant transmission and recovery rates are reasonable assumption at later stages. Computational results can be found in Figs. 13 and 14. We also notice that the number of infected people rises in summer, possibly due to more social contacts. This could eventually be regarded as the beginning of a ‘second wave’.

Figure 13

Computational results for our model of infected people in Germany for user-chosen , , for short-time simulation on (1 May 2020) and (1 September 2020) with real-world data as initial conditions

Figure 14

Computational results for our model of recovered people in Germany for user-chosen , , for short-time simulation on (1 May 2020) and (1 September 2020) with real-world data as initial conditions

For future investigations, it might be interesting to use more complex transmission rates or piecewise defined functions with switching points, see e.g. [49].

Computational short-time results for data from Iran

Real-world data from Iran and short-term computational results for COVID-19 data from Iran can be found in Figs. 15–23. Figure 17 again supports the assumption of an exponentially decaying transmission rate at the beginning of the spread of life-threatening disease. The computation results, depicted in Figs. 22 and 23, show qualitative agreement with the trends in real-world data. These results indicate that time-dependent transmission rates are a necessary addition to the classical SIR model. Alternatively, models with fractional derivatives could be considered [53, 54].

Figure 15

Unprocessed data for Iran with (1 March 2020) and (1 May 2020)

Figure 23

Computational results for our model with transmission and recovery rates from Figs. 20 and 21 of recovered people in Iran for short-time simulation on (1 March 2020) and (1 May 2020) with real-world data as initial conditions

Figure 17

Time-varying transmission rate from real-world data for Iran with (1 March 2020) and (1 May 2020)

Figure 22

Computational results for our model with transmission and recovery rates from Figs. 20 and 21 of infected people in Iran for short-time simulation on (1 March 2020) and (1 May 2020) with real-world data as initial conditions

Processed data for Iran with (1 March 2020) and (1 May 2020)

Time-varying recovery rate from real-world data for Iran with (1 March 2020) and (1 May 2020)

Time-varying basic reproduction number from real-world data for Iran with (1 March 2020) and (1 May 2020)

Figure 20

Time-varying transmission rates from real-world data for data from Iran with (1 March 2020) and (1 May 2020). The user-chosen parameters for our model are and

Figure 21

Recovery rates from real-world data for data from Iran with (1 March 2020) and (2 May 2020). The user-chosen recovery rate reads

Conclusion and outlook

We established certain properties such as well-posedness of the solution of our time-continuous SIR model in Sect. 2. Fortunately, we were able to transfer many properties of the time-continuous model to our time-discrete implicit SIR model in Sect. 3. These include unique solvability and monotonicity properties. In contrast to many other works mentioned in Sect. 1, we avoid an explicit forward model, but we could transform our implicit scheme to an easily solvable scheme. Thus, this makes our proposed scheme an attractive first prediction choice. In addition to that, we showed that our numerical scheme possesses an upper error bound.

Regarding our computational results, we see that our parametrization is an appropriate fit for first forecasts considering the first wave of a spreading virus. Since these transmission rates are monotonically decreasing, we, however, remark that we will need to use another parametrization if we want to model diseases with seasonal behavior [55]. Regarding our chosen examples, we get reasonable results. Additionally, we observe that our theoretical findings regarding monotonicity of recovered people from Theorem 10 are fulfilled in both examples. This stresses the attractiveness of our implicit solution scheme.

As depicted in Sect. 4, the inverse problem definitely needs further investigation. This is a topic of its own interest [56-59] since we need tools from different mathematical disciplines. As future research directions, extensions to further epidemiological forward models should be considered as we surely need more tools to predict the impact of upcoming epidemics. One can also consider delayed-differential or stochastic variants of our SIR model or modifications and extensions [23, 27] because, from a biological point of view, we often have to face integration of incubation times in epidemic models.