Stability analysis of a nonlocal SIHRDP epidemic model with memory effects.

The prediction and control of COVID-19 is critical for ending this pandemic. In this paper, a nonlocal SIHRDP (S-susceptible class, I-infective class (infected but not hospitalized), H-hospitalized class, R-recovered class, D-death class and P-isolated class) epidemic model with long memory is proposed to describe the multi-wave peaks for the spread of COVID-19. Based on the basic reproduction number R 0 , which is completely controlled by fractional order, the stability of the proposed system is studied. Furthermore, the numerical simulation is conducted to gauge the performance of the proposed model. The results on Hunan, China, reveal that R 0 < 1 suggests that the disease-free equilibrium point is globally asymptotically stable. Likewise, the situation of the multi-peak case in China is presented, and it is clear that the nonlocal epidemic system has a superior fitting effect than the classical model. Finally an adaptive impulsive vaccination is introduced based on the proposed system. Then employing the real data of France, India, the USA and Argentina, parameters identification and short-term forecasts are carried out to verify the effectiveness of the proposed model in describing the case of multiple peaks. Moreover, the implementation of vaccine control is expected once the hospitalized population exceeds 20 % of the total population. Numerical results of France, Indian, the USA and Argentina shed light on the varied effect of vaccine control in different countries. According to the vaccine control imposed on France, no obvious effect is observed even consider reducing human contact. As for India, although there will be a temporary increase in hospitalized admissions after execution of vaccination control, COVID-19 will eventually disappear. Results on the USA have seen most significant effect of vaccine control, the number of hospitalized individuals drops off and the disease is eventually eradicated. In contrast to the USA, vaccine control in Argentina has also been very effective, but COVID-19 cannot be completely eradicated.

Introduction

A dangerous infectious disease COVID-19 caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) was discovered at the end of 2019. Soon after, the number of confirmed cases continues to increase from day to day. As of January 1, 2021, a total of 83561252 cases have been diagnosed worldwide resulting in 1820668 deaths [1]. The rapid increase in infection cases indicates that the transmission capacity of COVID-19 is strictly stronger than that of MERS-CoV and SARS coronavirus [2, 3]. On March 11, 2020, COVID-19 has been declared a pandemic by World Health Organization. Meanwhile, multiple peaks of disease transmission have appeared around the world. For example, some countries are experiencing recurrent outbreaks. A better understanding and insight of the pandemic’s trend with the aid of mathematical model can provide a guideline for government and individuals to take appropriate actions to control the return of the epidemic.

The study of epidemiology, which deals with the dynamic evolution of diseases in a population, has attracted widespread attention in recent years [4-6]. The first model, introduced by Kermack and McKendrick, is SIR (susceptible-infected-recovered) model [7]. In this simplest SIR epidemic model, a basic assumption is that an individual’s time interval of staying in a compartment won’t have an influence on its transition out. This assumption is mathematically equivalent to a memoryless Markov process in which time is exponentially distributed at each interval. However, Angstmann et al. derived a fractional-order infected model, which may be adapted to any compartment model where the transition out of a compartment is dependent on the length of time since entering this compartment. Also, they pointed out that for some diseases with a chronic incubation period, such as human papiliomavirus (HPV), the distribution of infection time is a power-law tails [8]. Its behavior represents a much slower decay in time than an exponential decay, which is called nonlocal model with long memory. Thus, it is of much interest to consider nonlocal model to understand and control complex processes. Under this circumstance, many researchers established the nonlocal model, such as non-Markovian jump processes, Levy flights [9-11], fractional differential operators [12-14] and so on.

However, in the existing literature of nonlocal system, the kernel function is required to be integrable, or the kernel function is controlled by a decreasing exponential function, which results in system stability not being affected by delays. This suggests that the kernel function decays faster in time than the exponential distribution, see [15-22] and references therein. Zhang et al. investigated the dynamics of an SIR epidemic disease model with distributed times delay and the kernel function [15]. The dynamics of multi-group SEIR epidemic models with distributed and infinite delay are investigated by Shu et al. [16] and the kernel function . Furthermore, Muroya et al. successfully established a delayed multi-group SIS epidemic model with nonlinear incidence rates and patch structure, and they found the incidence time and the migration time did not alter the quality of the disease dynamics [21]. Hethcote et al. developed a epidemic model for which infection confers permanent immunity. The system was expressed by nonlinear Volterra integral equations of convolution type and the results showed that the distributed time also didn’t change the asymptotic behaviors of their model [22]. However, because COVID-19 has an incubation period of 7-14 days and is highly infectious during the incubation period, it is highly unreasonable to limit the kernel function. Meanwhile, little is known about COVID-19 as it is first transmitted in humans and different countries are going through this variety of multi-peak conditions, so policy formulation and action are strictly dependent on the history of disease transmission at different peak stages. Plus the precaution of people is also based on their knowledge on the primary course of COVID-19. In order to properly describe the multiple stages of this disease transmission, it is important to consider the kernel function of long-term memory. Few works have been devoted to studying the nonlocal epidemic model [23, 24]. Motivated by this, we focus on memory effects in this work, which means arbitrarily long history can be included.

Meanwhile, many countries have taken different measures to control the spread of COVID-19, such as sealing cities, wearing masks and keeping safe distances and so on. But as the epidemic continues to develop, the world is experiencing a resurgence of COVID-19. To scale back the possibilities of multiple transmission, the other suggested measures are the start of vaccination, the use of antiviral drugs, and the awareness about the risks of disease transmission, of which vaccination is considered to be one of the most effective way to control the spread of COVID-19. Thus, research into vaccines has become a paramount task. However, considering public opinion, diplomacy, economics, public health, and so on, an efficient vaccine allocation strategy appears to be equally compelling. Emanuel et al. offered a more ethically defensible and practical proposal for the fair distribution of COVID-19 vaccine: the Fair Priority Model [25]. Santiago et al. proposed a control strategy for SIRS epidemic model [26]. The model has been demonstrated to be stable for such vaccination strategy. Asier et al. presented a sliding mode controller for SEIR epidemic model [27]. Sultan et al. devoted to a stochastic model on the spread of COVID-19 and used the supermartingale approach to investigate a bound managing of which also leads to decrease of the number of infected individuals [28]. Rohith et al. proposed an appropriate threshold by a synthesis of government control variable, and a closed-loop control approach using sliding mode was adapted to control the COVID-19 pandemic [29].

Based on the above analysis, a nonlocal epidemic model with long memory is proposed for modeling the spread of COVID-19, where the infection of time conforms to a power-law tails distribution. Meanwhile based on the proposed model, an adaptive impulsive vaccination strategy with reference to the number of hospitalized populations is designed. The main contributions of this study are as follows:The rest of this paper is organized as follows: COVID-19’s SIHRDP nonlocal model with long memory is established in Sect. 2 and some preliminaries are also given in this section. Then, some dynamic behaviors of the proposed system are analyzed in Sect. 3. In Sect. 4, some dynamic behaviors of impulsive system are investigated. In Sect. 5, numerical simulations are provided to illustrate the theoretical results. Finally, the conclusions are given in Sect. 6.

The dynamic behaviors of the proposed system are analyzed by demonstrating the existence and uniqueness of the nonnegative solution, the global stability of the disease-free equilibrium point and the endemic equilibrium point.

Numerical results are presented to verify the theoretical results, in which the detection ratio is considered as a piecewise function based on the measures to strengthen nucleic acid detection.

Pulse vaccine control is investigated to prevent further spread of COVID-19. The basic reproduction number is obtained as a threshold condition to determine whether the disease will be stamped out.

Using the real data of France, India, the USA and Argentina, short-term forecasts of the number of hospitalized individuals at 7 days are given.

The adaptive pulse vaccine control is designed to update the vaccine rate in France, India, the USA and Argentina, when the hospitalized population exceeds of the total population.

Preliminaries and problem formulation

With the advent of the digital age, more and more attention has been paid to non-local models. This is one reason for the singular/anomalous behavior associated with complex systems and the increasing number of nonlocal phenomena. Before presenting the main result, some necessary preliminaries are introduced.

Preliminaries

This section begins with some definitions and results.

Definition 1

[30] A nonlocal linear operator is defined by the following form:where is the most general form of the kernel, as implied by the Schwartz kernel theorem [27].

Definition 2

[32] A Gamma function of is defined by:where .

Definition 3

[32] Let be an infinite interval on the real axis . The Riemann–Liouville fractional integrals of order is defined bywhere is the Gamma function.

Remark 1

When in Definition 1, one hasthat is, the Riemann–Liouville fractional integrals is a special nonlocal operator.

Lemma 1

[33] Suppose that a function is twice differentiable and that . Then, there are sequences and such thatand

Lemma 2

[32] Let F(s) be the Laplace transform of a function f(t), then the Laplace transform of the Riemann–Liouville fractional integral is given by

Lemma 3

Let and denote the Riemann–Liouville fractional integrals of a function f(t), then

Proof

The Laplace transform of the above terms is firstly considered as follows:Expand the above equation using a Taylor series expansion:Then the inverse Laplace variation of both sides can be obtained as follows:Therefore, it can be concluded from Eq. (2) as follows:

Lemma 4

[34] Consider the following impulsive differential equations:where , and . Then there exists a unique positive periodic solution of the above system:which is globally asymptotically stable, where

Problem formulation

For infectious diseases, information about disease history does not have the same impact all the time. Also the time it takes from incubation to infectivity may vary from person to person. Following the approach of [35] and taking the time of infection as a parameter, an SIHRDP epidemic model with distributed delay is considered as follows:The population density of a city is classified into S(t), I(t), H(t), D(t), R(t) and P(t), denoted susceptible, infective (infected but not hospitalized), hospitalized, dead, recovered and isolated class, respectively. The parameter is the inflow rate of susceptible individuals; is the transmission coefficient; is the recovery rate; is mortality rate; is the protection rate; is the disease-related mortality rate; is detection rate. Let K(s) denote the probability that an infected individual remains infectious and is the probability of surviving death from time s to time t, then the death survival can be written as . represents all the infected individuals I(s) before time t can survive to time t and remain infective. According to [36], assume that the probability of remaining in the infected class K(s) is the step-function given byfor finite . Then system (3) reduces to the following delay differential equation system:Moreover, the weak delay kernel function K(s) represents the impact of past memory on current dynamics and its form is as follows in [17]:where is the decline rate of the effects of past memories. Then it is easy to see that the dependence of system state on historical information is negligible as time goes by. However, when a disease spreads firstly between people, an individual’s prior experience or knowledge of the disease influences their subsequent response. An appropriate choice, in order to include long-term memory effects, could be a power-law function that shows a slow decay, so that the state of the previous state of the system also contributes to the evolution of the system. Motivated by this, consider the following power-law correlation function for K(s):in which and denotes the Gamma function. The choice of the coefficient and exponent allows us to rewrite system (3) to the form of Riemman–Liouville fractional integrals equations as follows:where is the Riemman–Liouvier fractional-order integral of . The decaying rate of the memory kernel (a time-correlation function) depends on . A smaller corresponds to a slower decay rate (long memory). Therefore, the length of memory is controlled by . As , the effect of memory is minimized: system (5) tends to be a memoryless system. A visual representation of the system structure can be seen in Fig. 1.

Fig. 1

A block diagram of the SIHRDP epidemic model for COVID-19

Model analysis

In this section, some dynamical behaviors of the proposed system (5) are analyzed. Here it can be found that the susceptible class S(t), the infected class I(t) and the hospitalization class H(t) are not affected by the recovered class R(t), the death class D(t) and the isolated class P(t) of system (5). Hence, we shall focus our attention on the following reduced system:

Nonnegative and boundedness

In this subsection, the existence, uniqueness and boundedness of a nonnegative solution for system (6) are considered.

Theorem 1

For any given initial condition , there exists a unique nonnegative bounded solution (S(t), I(t),  H(t) of system (6) for all .

According to the hypotheses stated by Miller [37], it is sufficient to ensure the existence and uniqueness of solutions of system (6). In addition, it can be concluded from the first equations in system (6) thatwhich implies that for . Next the nonnegativity of I(t) and H(t) will be studied. First assume that and claim both I(t) and H(t) are nonnegative for . If not, there must exist a finite time such that and for and either (i) and , or (ii) and . But when (i) holds, based on the second equation in system (6), one hasthen this implies thatwhich shows that (i) is impossible. Next assume that (ii) is true. As in (i), one has for . The third equations in system (6) give which is contradictory. Therefore, I(t) and H(t) are nonnegative for all . In the following, the boundedness of system (6) will be discussed. Let and . Let denote the total population. Add all equations in system (6), it can be yielded thatthen this implies that . It is obvious that S(t), I(t) and H(t) are bounded.

Local stability

In this subsection, the existence of positive equilibrium points and the locally stability analysis of system (6) will be discussed.

In general, the basic reproduction number is ‘the expected number of secondary cases from an infected individual in a fully susceptible population’ [41]. According to Watmough et al. [40], when the basic reproduction number , each individual releases on average less than one newly infected individual, and the infection cannot grow. Conversely, if , each infected individual produces more than one new infected person, and the disease will not die out. It is therefore important to describe the relationship between the basic reproduction number and the spread of infectious diseases. According to the basic reproduction number , the following results are obtained:

Theorem 2

If , system (6) has only a disease-free equilibrium point which is locally asymptotically stable. If , the disease-free equilibrium point is unstable.

It is easy to see that system (6) has only a disease-free equilibrium point . Then the stability of the disease-free equilibrium will be discussed. Linearizing system (6) at the disease-free equilibrium point gives the characteristic equation as follows:where . Then the characteristic equation can be reduced as follows:where . Since , it just needs to be shown that all roots of have negative real parts. Suppose where is a root of h(t), and rewrite asIt is concluded from thatwhere Q is calculated according to Lemma 3:Take the square of both sides of Eq. (8), one hasDenote and to be the left-hand and right-hand side of the above inequality, respectively. Then withandNotice thatandWhen , it is obviously to seewhich is a contradiction. Therefore, has no positive real root, and this implies that the disease-free equilibrium point is locally asymptotically stable for .

Assume . To show that the disease-free equilibrium point is unstable, it is only to show that admits a positive real root. Considering , it follows from the above equation thatLetthen one has by . Obviously, is an increasing function of z whereas is decreasing in terms of z. Thus, is an increasing function of z. Therefore, there must be a positive such that and this is a real positive root of h(z). Hence, the disease-free equilibrium point is unstable if .

Globally stability

In this subsection, the global stability of disease-free equilibrium point will be considered. Meanwhile, the existence and stability of the endemic equilibrium point will also be investigated.

Theorem 3

Consider system (6), if , the disease-free equilibrium point is globally asymptotically stable, else, the disease-free equilibrium point is unstable.

According to Theorem 2, it has been known that if , the disease-free equilibrium point is locally asymptotically stable. To this end, it is only to show it is globally attractive. Note that system (6) is bounded, that is S(t), I(t) and . Then denoteso , and exist with , and . By Lemma.1, there exists a sequence with as such thatThe second equation of system (6) can be rewritten asTake the limit of both sides of the above equations, one hasthat isIf and , the above equality yieldsOn the other hand, consider the first equation of system (6), one haswhich implies . This is a contradiction. So it can be concluded , that is . Hence, based on the third equation in system (6), one hasTake the limit of both side, one has , which implies . Similarly, it can be proved that . Therefore, the disease-free equilibrium point is globally asymptotically stable if . When , it can be concluded from Theorem. 2 that the disease-free equilibrium point is unstable.

Theorem 4

If , system (6) has a unique endemic equilibrium point which is globally asymptotically stable.

Take the equilibrium point to be such thatone hasIt is obvious that based on Lemma 3. Therefore, by solving the above equation, one has , and . As the population in each compartment cannot be negative, the endemic equilibrium point can only exist if . Therefore, when , system (6) has only a unique endemic equilibrium point. Define the Lyapunov functionswhereandThen the derivative of along the solutions of system (6) is given bySimilarly, based on the definition of , the derivative of along the solutions of system (6) is given byThen the derivative of along the solutions of system (6) is given byAdd Eqs. (11), (12) and (13), one haswhere for if and only if for . Therefore, it can be see that and only and if only and . Therefore, by LaSalle’s invariance principle, the endemic equilibrium point is globally asymptotically stable.

Impulsive vaccination

When COVID-19 was first spread, countries around the world took active measures to effectively contain the spread of the epidemic. But almost all countries are now experiencing a secondary spread of the disease. With unremitting efforts, the development of vaccines has been remarkably successful. But the question is how to properly delivery vaccines to control the spread of COVID-19. This section studies the dynamic behavior of system (6) under the regular impulsive vaccination.

Regular impulsive vaccination

The regular impulsive vaccination is characterized by a constant vaccination rate and a constant inter-vaccination time interval . For (), consider the following system:and for , consider the following impulsive vaccination:where is vaccination rate; is time intervals between two consecutive impulsive vaccinations. The dynamics of the auxiliary system (15) and (16) will be discussed in the following.

Theorem 5

The periodic solution of system (15) and system (16) is , where

When I(t) and H(t) are entirely absent permanently, i.e., for all . It can be reduced to the following system:By Lemma.4, one has that the periodic solution of system (15) is globally asymptotically stable.

Theorem 6

Let (S(t), I(t), H(t)) be any solution of system (15), then is locally asymptotically stable when , where

Consider the following system with the small amplitude perturbations:where U(t) satisfiesand is a unit matrix. Hence, the fundamental solution matrix isThe linearization of system (15) becomesThe stability of the periodic solution is determined by the eigenvalues ofwhich , and . Therefore, when ,that is , which implies is locally stable.

Remark 2

If (i.e., in the absence of impulsive vaccination), then the impulsive reproduction number becomes the standard, non-impulsive, reproduction number, i.e., , and implies that the stability at the disease-free equilibrium point when is globally asymptotically stable.

Remark 3

Obviously, the impulsive reproduction number R . As a result, there is no need to perform pulse vaccine control when the disease is guaranteed to be eradicated. However, when the disease is persistently spread, it can be eliminated by adding pulse control so as to reduce the basic reproductive number.

Adaptive impulsive vaccination

According to the analysis in the previous section, when the basic reproduction number is less than 1, the disease dies out. However, periodic pulse vaccine control will lead to a waste of vaccine or fail to achieve the control purpose, so more complex pulse vaccine strategy is studied in this subsection.

Like [40], adaptive pulse vaccine strategy updates the impulsive vaccination rate by using a rule in reference to the quantity of hospitality class. Based on Theorem 6, it can be known that the epidemic disease is guaranteed to be eradicated when . So identified parameters of the system (6) according to the real data of COVID-19, a fixed value for is chosen such that and , which means the disease does not die out and persists for a long time without vaccine control, but when all susceptible people are vaccinated, the disease is bound to be stamped out. Then a database of for , whereby taking Eqs. 17 and 18 into account.

The aim of this rule is to increase pulse vaccination rates, and if the number of hospitalizations exceeds a predetermined value, to reduce it. To this end, the law used for updating such as the vaccination rate at each vaccination time instant is given by:where is an auxiliary value given by:which is a predefined constant.

Simulation

Firstly, the validity of the theoretical results is verified based on real data. However, due to the different responses of different countries after the secondary transmission of COVID-19, for example, the disease was effectively contained after China experienced the second peak, but other countries such as France, the USA, India and Argentina, are facing further rapid spread of the disease. So this section will consider China and other countries separately, and control the spread of COVID-19 through an adaptive vaccine pulse strategy.

Data source

The data used in this paper are from the Johns Hopkins University Center for System Science and Engineering [1]. The Johns Hopkins University published data of accumulated and newly confirmed cases, recovered case and death case of COVID-19 from January 22, 2020. The data of COVID-19 cases of different city from January 23, 2020, to January 1, 2021, will be used.

The description of the stability results

In this subsection, on the basis of system (6), the fitting effect on the current hospital cases is discussed, and the accuracy of Theorem 3 is verified.

By analyzing the cure and mortality data of Hunan in China, the appropriate recovery rate function and mortality function is selected from Table. 1 according to [41]. Based on the real data of Hunan in China, parameter identification is given in Table 2. From Table 2, it is easy to see that the transmission rate is higher, which represents the highly infectious of COVID-19, but the protection rate indicates that individuals have a strong sense of self-prevention. Also, the fractional-order , which means that historical information has a very important influence on the current spread of disease. On the basis of Table 2, system (6) can effectively fit the confirmed cases, which can be seen in Fig. 2. Furthermore, it is found from Fig. 3 that the disease-free equilibrium point is globally asymptotically stable, which is consistent with Theorem 3.

Table 1

Recovery rate function and mortality function

λ(t)	a_11+e^a_2(t-a_3)	a_11+e^-a_2(t-a_3)	a_1+e^-a_2(t+a_3)
κ(t)	b_1e^b_2(t-b_3)+e^-b_2(t-b_3)	b_1e^-b_2(t-b_3)^2	b_1+e^-b_2(t+b_3)

Table 2

Parameter identification of Hunan

Parameters	Hunan(China)	Parameters	Hunan (China)
Λ	0.2431	ρ	0.7092
β	1.1525	α	0.5502
δ	0.1424	ϵ	0.3068
λ	a_11+e^-a_2(t-a_3)	[a_1,a_2,a_3]	[0.334, 0.0702, 32.6864]
κ	b_1e^b_2(t-b_3)+e^-b_2(t-b_3)	[b_1,b_2,b_3]	[0.0267, 0.09, 30]

Fig. 2

The number of hospitalization case

Fig. 3

On the basis of system (6), the simulation curves of susceptible individuals, infected individuals and hospitalized individuals (the reproduction number is )

The description of multi-wave

The Spread of COVID-19 in China

In order to further measure the long memory of the spread of COVID-19, the multi-peak in China will be considered in this subsection. In contrast to the long-term memory, the disease only relies on the current information in the transmission process and has no dependence on the historical information. In other words, the memoryless kernel function is considered [17]:Then system (5) becomes the local epidemic model as follows:COVID-19 was discovered in Wuhan, China at the end of 2019. But since the end of March, the second wave of epidemics occurred in Beijing, Heilongjiang, Guangdong and other cities in China. In addition, due to the implementation of the CT diagnosis method and the nucleic acid test method, the number of confirmed cases increased sharply and reached a peak in early February. Therefore, we assume that starting from 12 February, the diagnosis rate will reach and maintain its maximal value. However, since the end of June, the third wave of epidemics occurred in Beijing. According to the news by Beijing, the detection of nucleic acid reagents was increased from 17 June. Therefore, the diagnosis rate do not maintain the maximal value from 17 June, but subject to a new distribution. Like [42], the diagnosed period of two peaks and three peaks are described by the following piecewise function:andwhere is the initial diagnosis rate, is the diagnosis rate of individuals and are the decreasing rate of the diagnosed period from symptom onset to be diagnosed, is the maximum diagnosis rate of individuals, is 13 February, is 17 June. From Figs. 4, 5 and 6, it can be found that system (6) fits the data better than system (19). Also, it can be seen from Tables. 3 and 4 that in the case of Heilongjiang, Guangdong and Beijing, historical information on the spread of the disease is very important to the present, and the impact of historical information is different at different peak stages. For example, for Guangdong, the impact of the first phase is greatest, while for Beijing, the impact of the second stage is the most obvious. Here represents the fractional order of the peak of the wave in Tables. 3, 4.

Fig. 4

The number of confirmed and recovery cases in Heilongjiang with two peaks (Nonlocal (above), local (below)

Fig. 5

The number of confirmed and recovery cases in Guangdong with two peaks (nonlocal (above), local (below)

Fig. 6

The number of confirmed and recovery cases in Beijing with three peaks (nonlocal (above), local (below))

Table 3

The fractional-order in Heilongjiang and Guangdong

Fractional-order	ϵ_1	ϵ_2	ϵ_3
Heilongjiang	0.1547	0.2296	0.1978
Guangdong	3.2718×10^-7	0.4023	0.4582

Table 4

The fractional-order in Beijing

Fractional-order	ϵ_1	ϵ_2	ϵ_3	ϵ_4
Beijing	7.7659×10^-5	4.0589×10^-14	0.2902	0.9022

The spread of COVID-19 in France, India, USA and Argentina

As of 1 January, the global coronavirus death toll surpassed 1.8 million and the confirmed cases reached a total of 83561252. The USA, Argentina, India, France and so on are among the worst hit by COVID-19 across the world, where confirmed cases are still rising. In this subsection, based on the data from January 22, 2020, to January 1, 2021, parameter identification of the above four countries are carried out through system (6), respectively. Meanwhile, numerical simulations fit the real data, and short-term prediction from January 2 to January 8 are made to compare with the real reported data.

France reimposed a second blockade from October 30, 2020, in response to the rapid rebound of COVID-19. In the second outbreak, COVID-19 has spread at an alarming rate, surpassing even the most pessimistic predictions. Therefore, in order to better describe the spread of COVID-19 in France, the spread process is analyzed in two phases: January 22–October 29, 2020; October 30, 2020–January 1, 2021. Parameter identification is shown in Table.5. And it can be seen from Table.5 that regardless of the phase, the transmission rates in France are high, and the first phase is significantly higher than the second phase. Meanwhile, both the detection rates are very low in France, indicating that France still needs to make further efforts to effectively control the evaluation of the epidemic. But compared with the first stage, the self-protection rate of individuals in the second stage increased significantly, suggesting that citizens have fully realized the severity of COVID-19. Furthermore, the second stage relies more on history . And the mortality and recovery rates are about the same between the two stages. In addition, according to the results of parameter identification, data fitting results and prediction results are shown in Fig. 7 and Table. 6. At the same time, based on previous analysis, adaptive pulse control for the spread of COVID-19 in France is implemented. Set , that is, pulse control is implemented from 1 January when the hospital population exceeds of the total population. The results of the implementation of control are shown in Fig. 8, which shows that the implementation of vaccine control can effectively reduce the spread of the epidemic, but there will be an upward phase before the down phase. Whereas it also suggests that vaccine control will not end the epidemic any time soon. Furthermore, consider the importance of reducing human contact, such quarantine measures are divided into four scenarios, as shown in Table 7. The results of implementing these two measures are shown in Fig. 9, which shows that the implementation of quarantine and vaccine can lead to a rapid decline in the epidemic without any further upward phase, but the effect of adding quarantine is not obvious, no matter how strong the quarantine is. It is therefore suggested that restrictions on population movement can be relaxed at the same time as vaccine controls are in place. Moreover, if the outbreak is to be completely eradicated, France will have to take other measures to deal with the spread of COVID-19.

Table 5

Parameter identification of France

Parameters	France	Parameters	France
Λ_1	4.9035×10^5	ρ_1	4.0528×10^-4
β_1	0.9487	α_1	0.9158
δ_1	0.0484	ϵ_1	1.5329×10^-6
λ_1	a_1+e^-a_2(t+a_3)	[a_1,a_2,a_3]	[0.0007, 1, 54.6647]
κ_1	b_1+e^b_2(t+b_3)	[b_1,b_2,b_3]	[0.0003, 0.7831, 50.9326]
Λ_2	8.0876×10^6	ρ_2	0.1002
β_2	0.7724	α_2	0.7458
δ_2	0.0065	ϵ_2	4.4897×10^-8
λ_2	a_1+e^-a_2(t+a_3)	[a_1,a_2,a_3]	[0.0005, 1, 54.3428]
κ_2	b_1+e^b_2(t+b_3)	[b_1,b_2,b_3]	[0.0002, 0.7835, 51.0205]

Fig. 7

The number of isolated cases in France

Table 6

Predicted number of confirmed individuals of France

Date	Real data	Predicted data
1.2	2403680	2390075
1.3	2415841	2402164
1.4	2418734	2414224
1.5	2436782	2426249
1.6	2460267	2438236
1.7	2460267	2450183
1.8	2498130	2462087

Fig. 8

The number of confirmed cases in France with vaccination and the change of vaccination rate

Table 7

The public measures are taken from the France COVID-19 alert system

Level	β (%)	Public measures parameters
1(Do nothing)	100	Vaccination and No public measures
2 (Prevent)	50	Vaccination, School and workplaces open
3 (Restrict)	25	Vaccination and Learn at home
4 (Lock-down)	12.5	Vaccination and Instructed to stay at home

Fig. 9

The number of confirmed cases in France with vaccination and quarantine (, , , respectively)

Meanwhile, as of 1 January, India had the second largest number of confirmed cases in the world with a cumulative total of 10256701. Based on real data from January 22, 2020, to January 1, 2021, parameter identification is shown in Table.8. As can be seen from Table.8, the sense of self-protection of Indian citizens is stronger than that of France and the transmission rate is strictly lower than that of France. Moreover, the transmission of COVID-19 in India is strictly dependent on historical information and cannot be considered as a memoryless process. The corresponding fitting results and forecast results are shown in Fig.10 and Table.9. Similarly, based on previous analysis, adaptive pulse control for the spread of the outbreak in India from January 1, 2021, is implemented. Set , that is, pulse control is implemented when the inpatient population exceeds of the total population. As can be seen from Fig. 11, pulse control also lead to an increase in the number of diagnoses, but the increase is slight. And, unlike the cases of France, vaccine control in India can quickly reduce the spread of COVID-19 and thus achieve eradication. What’s more, adaptive pulse control ultimately results in the proportion of vaccine 0.203 that makes the disease disappear.

Table 8

Parameter identification of India

Parameters	India	Parameters	India
Λ	6.8481×10^7	ρ	0.2251
β	0.3968	α	0.0979
δ	0.0061	ϵ	6.6656×10^-6
λ	a_11+e^-a_2(t-a_3)	[a_1,a_2,a_3]	[0.0912, 0.5222, 0]
κ	b_1e^b_2(t-b_3)+e^-b_2(t-b_3)	[b_1,b_2,b_3]	[0.0024, 0.0080, 0.0108]

Fig. 10

The number of confirmed cases in India

Table 9

Predicted number of confirmed individuals of India

Date	Real data	Predicted data
1.2	247220	246728
1.3	243953	243086
1.4	231236	239492
1.5	227546	235945
1.6	228083	232445
1.7	225449	228991
1.8	225449	225583

Fig. 11

The number of confirmed cases in India with vaccination and the change of vaccination rate

Furthermore, the USA is now the worst affected country by COVID-19 in the world, with 220,000 new confirmed cases in a single day. It becomes necessary to study how the epidemic will spread in the USA and predict the number of cases to be confirmed, as well as to come up with appropriate solutions. Moreover, the real data are given by [1], which stopped collecting recovery data from December 14 in the USA. And the number of confirmed individuals in the USA has increased dramatically since October 31, 2020. This may be due to changes in government policies and individual behavior, leading to a surge in the number of infected people. In order to respond effectively to the COVID-19 outbreak in the USA, the spread of the outbreak in the USA is analyzed in two phases: January 22, 2020–October 31, 2020; November 1, 2020–December 13, 2020. Therefore, based on the data from January 22, 2020, to December 13, 2020, the code like [41] is used for parameter identification, and the results are shown in Table. 10. As can be seen from Table.10, the transmission rate of COVID-19 in the USA is strictly lower than that in France and India, but the lower self-protection rate of American citizens, coupled with insufficient detection efforts, making the disease spread wantonly in the USA. In addition, the spread of the disease in the USA is less dependent on historical information, which may be due to different measures taken in the USA at different times. Furthermore, Table. 11 presents the short-term prediction from December 14 to December 20. Compared with the real reported data, it can be found from Fig. 12 that system (6) has a good fitting and forecasting effect. Similarly, based on previous analysis, we implemented adaptive pulse control for the spread of COVID-19 in the USA from December 14, 2020. As can be seen from Fig. 13, different from in France and India, vaccine control in the USA can quickly reduce the spread of the epidemic and thus achieve eradication. What’s more, adaptive pulse control ultimately results in the proportion of vaccine 0.243 that makes the disease disappear.

Table 10

Parameter identification of the USA

Parameters	US	Parameters	US
Λ	9.8523×10^6	ρ	0.0249
β	0.783	α	0.7725
δ	0.0236	ϵ	8.5557×10^-6
λ	a_11+e^-a_2(t+a_3)	[a_1,a_2,a_3]	[0.008, 0.0497, 137.6922]
κ	b_1+e^-b_2(t+b_3)	[b_1,b_2,b_3]	[0.0003,0.0438,163.1724]
Λ	5.5672×10^6	ρ	0.0088
β	0.1841	α	0.3773
δ	0.017	ϵ	0.7311
λ	a_11+e^-a_2(t+a_3)	[a_1,a_2,a_3]	[0.0533, 0.0082, 295.1663]
κ	b_1+e^-b_2(t+b_3)	[b_1,b_2,b_3]	[0.0002, 0.4180, 17.7102]

Table 11

Predicted number of confirmed individuals of the USA

Date	Real data	Predicted data
12.14	16279098	15373133
12.15	16474503	15516735
12.16	16718162	15701502
12.17	16952589	15849467
12.18	17199987	15891954
12.19	17393546	15981480
12.20	17581534	16073154

Fig. 12

The number of confirmed cases in the USA

Fig. 13

The number of confirmed cases in the USA with vaccination and the change of vaccination rate

Finally, to further verify the effect of fractional order on disease development, the transmission of COVID-19 in Argentina is considered based on data from January 22, 2020, to January 1, 2021. Using the method similar to the other three countries, parameter identification is obtained in Table.12, and the transmission of COVID-19 in Argentina is strictly dependent on historical information and the effect of memory length is lower than in the other three countries. The corresponding fitting results and forecast results are shown in Fig.14 and Table.13. Similarly, based on previous analysis, adaptive pulse control for the spread of the outbreak in Argentina from January 1, 2021, is implemented. As can be seen from Fig. 15, pulse control reduces the number of confirmed cases immediately, but the outbreak can not be eradicated. At the same time, the proportion of vaccine delivered is not constant, which fluctuates up and down.

Table 12

Parameter identification of Argentina

Parameters	Argentina	Parameters	Argentina
Λ	1.9389×10^6	ρ	0.039
β	0.3583	α	0.3049
δ	0.0039	ϵ	0.0939
λ	a_1+e^-a_2(t+a_3)	[a_1,a_2,a_3]	[0.0161, 0.0127, 207.8303]
κ	b_1e^b_2(t-b_3)+e^-b_2(t-b_3)	[b_1,b_2,b_3]	[0.0045, 0.0214, 2.4528]

Fig. 14

The number of confirmed cases in Argentina

Table 13

Predicted number of confirmed individuals of Argentina

Date	Real data	Predicted data
1.2	144367	155101
1.3	144276	156772
1.4	147223	158480
1.5	152329	160221
1.6	158147	161992
1.7	161090	1163789
1.8	164283	165610

Fig. 15

The number of confirmed cases in Argentina with vaccination and the change of vaccination rate

Conclusion and discussion

Considering the multi-peak situation of epidemic transmission, the nonlocal SIHRDP epidemic model with memory effect is investigated in this paper. Firstly, the existence and uniqueness of the nonnegative solution for system (6) are established. Then the global stability of the disease-free equilibrium point is obtained when the basic reproduction number , and when , the global stability of the endemic equilibrium point is also discussed. Meanwhile, with the real data from January 22, 2020, to January 1, 2021, in Hunan, parameter identification is carried out and the basic reproduction number , which represents COVID-19 would not spread in Hunan if the current policy remained unchanged. In addition, the fractional-order indicates that the spread of the disease in Hunan is strictly dependent on historical information. Besides, in order to verify the validity of system (6), based on the real data of Beijing, Guangdong and Heilongjiang from January 22, 2020, to January 1, 2021, multiple peaks in China are considered. Since the Chinese government increased nucleic acid testing, the diagnosis rate is considered as a piecewise function. After fitting the real data through system (6) and system (19), the results reveal that the transmission of disease in these three provinces has shown the memory effect.

At the same time, with the success of the vaccine development, how to distribute the vaccine reasonably becomes very important. In this paper, adaptive pulse vaccine control is considered based on system (6). Firstly, the basic reproduction number under the pulse condition is obtained, and it is proved that if , the disease will die out. Then, employing the basic reproduction number, an adaptive pulse control is designed to control the outbreak as much as possible. Meanwhile, parameter identification is carried out with real data in France from January 22, 2020, to January 1, 2021. Also it can be seen that the transmission rate is very high; however, the self-protection awareness of citizens is very slack, and the testing is not sufficiently enforced. Moreover, no matter in which stage of epidemic transmission, there is a strict dependence on historical information. Moreover, based on the parameter identification results, the short-term prediction of COVID-19 in France is made, and it is found that system (6) not only has a good fitting effect, but also has a good prediction effect. Then we presumed that when the number of hospitalizations exceeds 20 of the total population, vaccine control will be implemented. The numerical results show that the number of confirmed cases in France will experience a brief rise and then a decline, but COVID-19 will not disappear completely. Then, by setting up four quarantine measures, it can be found that reducing exposure under vaccine control conditions not doesn’t make a radical difference. Therefore, based on the above analysis, quarantine measures can be relaxed under vaccine control; moreover, other effective measures must be taken in France to eliminate COVID-19. Compared with France, India has a strong sense of self-protection among its citizens, and the rate of disease transmission in India is strictly lower that of France. Similarly, applying adaptive pulse control to India under the same conditions, it can be found that the number of confirmed also experiences a brief rise, but eventually falls to zero. Varied from the cases in France and India, the number of confirmed cases in the USA drops to zero as soon as vaccine control is implemented. But using the adaptive pulse vaccine to control outbreaks in Argentina, like those in the USA, would immediately decline, but not ultimately eradicate. In this paper, the methods used are the least square method and the first difference method, and the research on the optimization of the algorithm [43] will be our future work.

This paper presents the relationship between memory and disease transmission, and the numerical simulation confirms the derived theory by harnessing real data of COVID-19. Although this study has a high accuracy in predicting the epidemic situation and formulating reasonable measures, it also raises several questions for further research: how to set up reasonable measures in France to completely eliminate the epidemic situation; what to do when the effect and quantity of vaccine is overestimated; how to set the time arrangement of vaccine delivery more reasonable and so on.