ABSTRACT
Background The threats of the COVID-19 pandemic require the mobilization of scientists, including mathematicians. To understand how the number of cases increases versus time, various models based on direct observations of a random number of new cases and differential equations can be used. Complex mathematical models contain many unknown parameters, the values of which must be determined using a limited number of observations of the disease over time. Even long-term monitoring of the epidemic may not provide reliable estimates of its parameters due to the constant change of testing conditions, isolation of infected and quarantine. Therefore, simpler approaches should also be used, for example, some smoothing of the dependence of the number of cases on time and the known SIR (susceptible-infected-removed) model. These approaches allowed to detect the waves of pandemic in different countries and regions and to make adequate predictions of the duration, hidden periods, reproduction numbers, and final sizes of its waves. In particular, seven waves of the COVID-19 pandemic in Ukraine were investigated.
Objective We will detect new epidemic waves in Ukraine that occurred after September 1, 2020 and estimate the epidemic characteristics with the use of generalized SIR model. Some predictions of the epidemic dynamics will be presented.
Methods In this study we use the smoothing method for the dependence of the number of cases on time; the generalized SIR model for the dynamics of any epidemic wave, the exact solution of the linear differential equations and statistical approach developed before.
Results Seventh and eights epidemic waves in Ukraine were detected and the reasons of their appearance were discussed. The optimal values of the SIR model parameters were calculated. The prediction for the COVID-19 epidemic dynamics in Ukraine is not very optimistic: new cases will not stop appearing until June 2021. Only mass vaccination and social distancing can change this trend.
Conclusions New waves of COVID-19 pandemic can be detected, calculated and predicted with the use of rather simple mathematical simulations. The expected long duration of the pandemic forces us to be careful and in solidarity.The government and all Ukrainians must strictly adhere to quarantine measures in order to avoid fatal consequences.
Introduction
The studies of the COVID-19 pandemic dynamics in Ukraine are presented in [1-14] and summarized in the book [15]. Different simulation and comparison methods together with the official WHO data sets, [16] were used. In particular the classical SIR model [17-19], connecting the number of susceptible S, infected and spreading the infection I and removed R persons, was applied in [2-5, 8-10] to simulate the first pandemic wave in Ukraine. The unknown parameters of this model were be estimated with the use of the cumulative number of cases V=I+R and the statistics-based method of parameter identification developed in [20, 21].
The weakening of quarantine restrictions, changes in the social behavior and probably in the coronavirus activity causes change in SIR characteristics and the epidemic dynamics. To detect these changes, a simple method of numerical differentiations of accumulated number of cases was proposed in [11, 15]. To simulate these new pandemic waves, the SIR model was generalized in [12, 15]. In [12, 14, 15] the results of simulation of the first six epidemic waves in Ukraine are presented with the use of a procedure for sequentially determining the parameters of the model for each epidemic wave, starting with the first one.
This method requires considerable effort and time. The book [15] introduced a new algorithm for determining the optimal parameter values for a particular epidemic wave without calculating the dynamics of previous waves and presented calculations for the seventh epidemic wave in Ukraine. In this paper, we will analyze the dynamics of the epidemic in Ukraine in the period from September 1 to December 20, 2020, calculate the parameters of the eighth wave and make some predictions.
Materials and Methods
Data
The official information regarding the accumulated numbers of confirmed COVID-19 cases Vj in Ukraine according to the official sources [22, 23] is shown in Table 1. Unfortunately, WHO stopped to present the daily number of new cases in August 2020. The corresponding moments of time tj (measured in days, zero point is January 20, 2020) are also shown in this table. To calculate the SIR characteristics for the seventh epidemic wave in Ukraine, the data set for the period October 1-14 was used in [15]. We will use the period November 21 – December 4, 2020 to calculate the characteristics of the eighth epidemic wave. Other values presented in Table 1 were used only for comparisons and verifications of calculations.
It must be noted that the data presented in Table 1 does not show all the COVID-19 cases in Ukraine. Many infected persons are not identified, since they have no symptoms. Many people know that they are ill, since they have similar symptoms as other members of families, but avoid to make tests. Unfortunately, one laboratory confirmed case can correspond to several other cases which are not confirmed and displayed in the official statistics. This fact reduces the accuracy of mathematical simulations.
Detection of epidemic waves
To control the changes of epidemic parameters, we can use daily numbers of new cases and their derivatives. Since these values are random, we need some smoothing. For example, we can use the smoothed daily number of accumulated cases proposed in [11, 12, 14, 15]:
The first and second derivatives can be estimated with the use of following formulas:
Generalized SIR model
The classical SIR model for an infectious disease [17-19] was generalized in [12, 15] in order to simulate different epidemic waves. We suppose that the SIR model parameters are constant for every epidemic wave, i.e. for the time periods: . Than for every wave we can use the equations, similar to [17-19]:
Here S is the number of susceptible persons (who are sensitive to the pathogen and not protected); I is the number of infected persons (who are sick and spread the infection; please don’t confuse with the number of still ill persons, so known active cases) and R is the number of removed persons (who no longer spread the infection; this number is the sum of isolated, recovered, dead, and infected people who left the region). Parameters αi and ρi are supposed to be constant for every epidemic wave.
Parameters αi show how quick the susceptible persons become infected (see (4)). Large values of this parameter correspond to severe epidemics with many victims. These parameters accumulate many characteristics. First they shows how strong (virulent) is the pathogen and what is the way of its spreading. Parameters αi accumulate also the frequency of contacts and the way of contacting. In order to decrease the values of αi, we have to minimize the number of our contacts and change our contacting habits. For example, we have to avoid the public places and use masks there, minimize or cancel traveling. We have to change our contact habits: to avoid handshakes and kisses. First, all these simple things are very useful to protect yourself. In addition, if most people follow these recommendations, we have chance to diminish the values of parameters αi and reduce the negative effects of the pandemic.
The parameters ρi characterize the patient removal rates, since eq. (6) demonstrates the increase rate of R. The inverse values 1/ ρi are the estimations for time of spreading infection τi during i-th epidemic wave. So, we are interested in increasing the values of parameters ρi and decreasing 1/ ρi. People and public authorities should work on this and organize immediate isolation of suspicious cases.
Since the derivative d (S + I + R) / dt is equal to zero (it follows from summarizing Eqs. (4)-(6)), the sum must be constant for every wave and is not the volume of population.
To determine the initial conditions for the set of equations (4)–(6), let us suppose that at the beginning of every epidemic wave :
It follows from (4) and (5) that
Integration of (9) with the initial conditions (8) yields:
It follows from (9) that function I has a maximum at S =vi and tends to zero at infinity. The corresponding number of susceptible persons at infinity Si∞ > 0 can be calculated from a non-linear equation
Formula (11) follows from (10) at I=0.
In [12, 15] the set of differential equations (4)-(7) was solved by introducing the function corresponding to the number of victims or the cumulative confirmed number of cases. For many epidemics (including the COVID-9 pandemic) we cannot observe dependencies S(t), I (t) and R(t) but observations of the accumulated number of cases Vj corresponding to the moments of time tj provide information for direct assessments of the dependence V (t).
It follows from (5) and (6) that:
Eqs. (7), (10) and (13) yield:
Integration of (14) provides an analytical solution for the set of equations (4)–(6):
Thus, for every set of parameters Ni, Ii, Ri, vi, αi and a fixed value of V, integral (16) can be calculated and the corresponding moment of time can be determined from (15). Then functions I(t) and R(t) can be easily calculated with the use of formulas (10) and:
The final numbers of victims (final accumulated number of cases corresponding to the i-th epidemic wave) can be calculated from:
To estimate the final day of the i-th epidemic wave, we can use the condition: which means that at t > tif less than one person still spreads the infection.
Parameter identification procedure
In the case of a new epidemic, the values of its parameters are unknown and must be identified with the use of limited data sets. For the first wave of an epidemic starting with one infected person, the number of unknown parameters is only four, since I1 = 1 and R1 = 0. The corresponding statistical approach was proposed in [20, 21] and used [2-5, 8-10] to simulate the first COVID-19 pandemic wave in Ukraine.
For the next epidemic waves (i > 1), the moments of time corresponding to their beginning are known. Therefore the exact solution (15)-(16) depend only on five parameters - Ni, Ii, Ri, vi, αi. Then the registered number of victims Vj corresponding to the moments of time tj can be used in eq. (16) in order to calculate for every fixed values of Ni,vi, Ii, Ri and then to check how the registered points fit the straight line (15).
Eq. (15) can be rewritten as follows:
Assuming we can estimate the values of parameters γ and β, by treating the values and corresponding time moments tj as random variables Then we can use the observations of the accumulated number of cases and the linear regression in order to calculate the coefficients and of the regression line using the standard formulas from, e.g., [24]. Values and can be treated as statistics-based estimations of parameters γ and β from relationships (21).
The reliability of the method can be checked by calculating the correlation coefficients ri (see e.g., [24]) for every epidemic wave checking how close its value is to unity. We can use also the F-test for the null hypothesis that says that the proposed linear relationship (20) fits the data set. The experimental values of the Fisher function can be calculated for every epidemic wave with the use of the formula: where ni is the number of observations for the i-th epidemic wave, m = 2 is the number of parameters in the regression equation. The corresponding experimental value Fi has to be compared with the critical value FC (k1, k2) of the Fisher function at a desired significance or confidence level (k1 = m −1, k2 = ni − m). When the values ni and m are fixed, the maximum of the Fisher function coincides with the maximum of the correlation coefficient. Therefore, to find the optimal values of parameters Ni, vi, Ii, Ri, we have to find the maximum of the correlation coefficient for the linear dependence (20). To compare the reliability of different predictions (with different values of ni) it is useful to use the ratio Fi / FC (1, ni − 2) at fixed significance level. We will use the level 0.001; corresponding values of FC (1, ni − 2) can be taken from [25]. The most reliable prediction yields the highest Fi / FC (1, ni − 2) ratio.
The exact solution (15)-(16) allows avoiding numerical solutions of differential equations (4)-(6) and significantly reduce the time spent on calculations. In the case of sequential calculation of epidemic waves i = 1,2,3 …, it is possible to avoid determining the four optimal unknown parameters Ni, vi, Ii, Ri, thereby reducing the amount of calculations and difficulties in isolation a maximum of the correlation coefficient. For parameters Ii, Ri it is possible to use the numbers of I and R calculated for the previous wave of epidemic at the moment of time when the following wave began. Then we need to calculate values , linear regression coefficients (22), correlation coefficient ri, Fi / FC (1, n − 2) and to isolate the values of parameters Ni and v I corresponding to the maximum of ri. Knowing the optimal values of five parameters Ni, Ii, Ri, vi, αi, the SIR curves and other characteristics of the corresponding epidemic wave can be calculated with the use of formulas (10)-(17). This approach has been successfully used in [12, 14, 15]. In particular, six waves of the Covid-19 epidemic in Ukraine and four pandemic waves in the world were calculated.
Segmentation of epidemic waves and their sequential SIR simulations need a lot of efforts. To avoid this, a new method of obtaining the optimal values of SIR parameters was proposed in [15]. First of all we can use the relationship which follows from (12). To estimate the value Vi, we can use the smoothed accumulated number of cases (e.g., formula (1)). Then where i corresponds to the moment of time . To obtain one more relationship, let us use (7) and (13)
To estimate the average number of new cases dV/dt at the moment of time , we can use (2). Thus we have only two independent parameters Ni and vi. To calculate the value of parameter αi, some iterations can be used (see details in [15]).
Results and discussion
The COVID-19 pandemic characteristics for Ukraine in autumn 2020 are shown in Fig. 1. Differentiation of the smoothed number of accumulated cases (eq. (1), line) with the use of formulas (2) (“triangles”) and (3) (stars) allow us to detect the changes in epidemic dynamics. It can be seen that after November 15 the average daily number of new cases (“triangles”) started to decrease. Similar short periods of the epidemic stabilization occurred in May, June and August, 2020 (see [11, 12, 14, 15]). Let us hope that the Christmas and New Year celebrations will not significantly worsen the existing positive dynamics.
The values of the second derivative (3) allows to detect the changes in the epidemic characteristics and to separate its different waves. The jump in d2V/dt2 values corresponding to September 11 (see “stars” in Fig. 1) can be explained by the beginning of classes in schools and universities (on September 1). Children and young people are often asymptomatic carriers of the infection and bring it to their families. For example, employees of two kindergartens and two schools in the Ukrainian city of Chmelnytskii were tested for antibodies to COVID-19, [26]. In total 292 people work in the surveyed institutions. Some of the staff had already fallen ill with COVID-19 or were hospitalized. Therefore, they were not tested accordingly. Of the 241 educators tested, antibodies were detected in 148, or 61.4%. These results indicate the important role of children in the spread of COVID-19 infection and the fact that in Ukraine those people who have become ill and have antibodies to coronavirus infection, obviously, are much more than in the official recovery statistics presented in Table 1.
The severe jumps in d2V/dt2 values occurred also in October and November, 2020 (see “stars” in Fig. 1). Probably, this is due to the local elections and a presidential poll, which were held throughout Ukraine on October 25, 2020 and involved hundreds of thousands of people to campaign and work in election commissions (their number was about 30 thousand). This obviously increased the number of contacts and the likelihood of additional infections. The corresponding seventh epidemic wave in Ukraine was considered in [15] (the results for the first six waves can be found in [12, 14, 15]). The optimal values of SIR parameters; predictions of the epidemic wave duration tif and its final size Vi∞ are shown in Table 2. The characteristics of the eights wave (occurred after November 21, 2020) were calculated and presented in Table 2. It can be seen that the predictions for Ukraine are very pessimistic. The number of cases will exceed 1.07 million and new cases will appear even in June 2021.
Fig. 2 illustrates the COVID-19 pandemic waves in Ukraine calculated with the use of the optimal values of SIR parameters presented in Table 2 and formulas (10)-(17). Solid lines correspond to the numbers of victims V=I+R, dashed lines show the numbers of persons spreading the infection I multiplied by 10, dotted lines illustrate the theoretical estimations of the daily numbers of new cases dV/dt calculated with the use of eq. (14) and multiplied by 10. The red, black and blue colors correspond to sixth, seventh and eighth epidemic waves respectively. The accumulated number of cases Vj taken for calculations are shown by circles; “triangles” and “stars” correspond the cases taken only for comparisons and verifications of the predictions.
Fig. 2 illustrates rather good accuracy of SIR simulations. In particular, the Vj values (“triangles”) started to deviate from the solid red line only after September 16. This fact can be explained by the beginning classes in schools and universities (as mentioned above). This curve was calculated with the use of data set from the period August 9-22. The increased number of contacts in schools and universities caused rapid increase of the number of cases. This increase in the number of diseases intensified in October and November 2020 through elections and a presidential poll. The very irregular nature of the epidemic dynamics (particularly the large values of the second derivative (see “stars” in Fig. 1)) led to the fact that V=I+R curve for the seventh wave very quickly began to deviate from the recorded number of cases Vj (compare solid black line and “stars” in Fig.2).
The absence of sharp jumps of the second derivative after November 21, 2020 (see “stars” in Fig. 1), allowed us to predict quite accurately the further epidemic dynamics (compare solid blue line and “stars” in Fig.2). Blue dashed line in Fig.2 show that the number of persons spreading the infection diminished in December 2020. Let us hope that the Christmas and New Year celebrations will not significantly worsen this positive dynamics.
Conclusions
Constant changes in the Covid-19 pandemic conditions (i.e., in the peculiarities of quarantine and its violation, in situations with testing and isolation of patients, in coronavirus activity due to its mutations etc.) lead to new pandemic waves and cause the changes in the values of parameters of the mathematical models. To identify these changes, a simple method was proposed based on the numerical differentiation of smoothed dependences of the accumulated number of cases. The results show that smoothing dependence of the accumulated number of cases and its differentiation can provide fairly accurate and useful information about the course of the epidemic, identify important changes in its dynamics and provide timely recommendations for quarantine measures or control of social distancing.
To simulate different pandemic waves (periods with more or less constant values of its dynamics parameters) a new generalized SIR model and its exact solution have been proposed. Procedure of its parameters identification were developed and successfully applied to calculate the characteristics of several pandemic waves in Ukraine.
The pandemic duration projections are unfortunately not optimistic. We hope that effective quarantine measures and successful vaccination will be able to reverse this trend. Very long duration of the pandemic requires correction of our behavior, we can not live as before it occurred. Decreased feelings of insecurity and non-compliance with social distancing may further increase the pandemic duration and the number of the coronavirus victims. Total closure of settlements or regions can be recommended only in the event of a sharp increase in the number of cases. There are many things that can be done without loss to the economy and our daily lives:
Minimize the number of contacts and trips, not visit crowded places. Work and study remotely where possible.
Refrain from shaking hands and kisses during meetings. Use masks in transport and crowded areas.
If you (or others) have any suspicious symptoms, do your best to avoid the spread of the infection.
Data Availability
All data are in the text
Acknowledgements
I would like to express my sincere thanks to Professors Dirk Langemann (Techniche Universitaet Braunschweig) and Juergen Prestin (Universitaet zu Luebeck) for their support in developing the used optimization approach. I would like to thank also Professors Alberto Redaelli, Giuseppe Passoni and Gianfranco Fiore (Politecnico di Milano), S. Pereverzyev (RICAM, Linz, Ausria) for involving me in very interesting biomedical investigations in frames of EU-financed Horizon-2020 projects EUMLS (Grant agreement PIRSES-GA-2011-295164-EUMLS) and AMMODIT (Grant Number MSCA-RISE 645672).
I would also like to thank Academician Viktor Grinchenko, Professors Volodymyr Tymofeyev, Alexander Galkin, Pavlo Maslianko; my friends Liudmyla Trotsenko, Maia Troian, Nina Basiuk, Anna Voronka, Daria Cusitcaia, Damien Berezenko, Anatolii Podkur, Volodymyr Borysenko and Oleksii Rodionov for their support and help in collecting and processing data.