Stability Analysis of an Extended SEIR COVID-19 Fractional Model with Vaccination Efficiency

This work is aimed at presenting a new numerical scheme for COVID-19 epidemic model based on Atangana-Baleanu fractional order derivative in Caputo sense (ABC) to investigate the vaccine efficiency. Our construction of the model is based on the classical SEIR, four compartmental models with an additional compartment V of vaccinated people extending it SEIRV model, for the transmission as well as an effort to cure this infectious disease. The point of disease-free equilibrium is calculated, and the stability analysis of the equilibrium point using the reproduction number is performed. The endemic equilibrium's existence and uniqueness are investigated. For the solution of the nonlinear system presented in the model at different fractional orders, a new numerical scheme based on modified Simpson's 1/3 method is developed. Convergence and stability of the numerical scheme are thoroughly analyzed. We attempted to develop an epidemiological model presenting the COVID-19 dynamics in Italy. The proposed model's dynamics are graphically interpreted to observe the effect of vaccination by altering the vaccination rate.


Introduction
Emerging in December 2019, COVID-19 became the biggest global threat. The whole world faced a critical situation, and a panic was created throughout the globe. The World Health Organization (WHO) declared it as an international health emergency. Not much was known about this virus as it was a new virus that has not been identified earlier in humans. So, relying on the observations and symptoms shown by the people who caught the virus, strategies of isolation, lockdown, and social distancing were adopted to prevent the disease from spreading. With the passage of time, the virus becomes more threatening because of its new variants. To investigate the dynamics and transmission of the virus, researchers utilized the concept of mathematical modeling [1][2][3][4]. These models actually do not provide the cure for the infectious diseases but can contribute in predicting the dynamics and behaviour of these infections through simulations.
Researchers rely more on fractional models because of their nonlocal behaviour, hereditary properties, and memory effects rather than integer order models. Besides, a number of experimental facts exhibit that the natural dynamics go along with fractional calculus. Ahmad et al. [5] utilized the modified Euler's method to study the dynamics of the model for COVID-19 based on fractional differential equations and simulate the results by using the available data of initial days of the spread in Wuhan City. Zakary et al. [6] investigated the situation of the spread of the Corona Virus under the effect of two classes, quarantined and other who did not respect the quarantine strategy using the mathematical model, and an optimal control strategy to reduce the infections in Morocco. Ahmad et al. [7] presented a brief analysis on the WHO reported data for Pakistan and investigated it using the SEIR fractional model. Askar et al. [8] considered the SITR fractional model to forecast the transmission of deadly virus along with imposed lockdown in India. Bushnaq et al. [9] investigated the effect of physical distance to control the spread of the virus through fractional order model. Abdulwasaa et al. [10] proposed a fractal fractional model for forecasting deaths and new cases of COVID-19 outbreak in India. Baba and Rihan [11] used the Caputo-Fabrizio fractional order model to analyze the dynamics of COVID-19 variants. A fractional order delay differential model for COVID-19 infection was provided by Rihan and Gandhi [12] to understand what triggers the severity of symptoms and disease of polluted lung and respiratory system. A stochastic SIAQR epidemic model was proposed by Rihan and Alsakaji [13] to understand the dynamics of COVID-19. Furati et al. [14] tested the data taken from China, Saudi Arabia, Brazil, and Italy by proposing a numerical scheme using fractional order model performing the two-step generalized exponential timedifferencing method. Alqahtani [15] numerically analyzed and investigated the stability of the SIR epidemic model considering the health system. Recently another fractional order derivative with nonsingular kernel named as Generalized Hattaf Fractional (GHF) derivative [16] is introduced. We can see the use of this derivative for epidemiological models in recent literature. By creating a new numerical method, Hattaf [17] investigated the qualitative properties of solutions to fractional differential equations with the new GHFs derivative, including stability, asymptotic stability, and Mittag-Leffler stability and later applied the obtained analytical and numerical results to a biological nonlinear system derived from epidemiology. Cheneke et al. [18] developed and examined a coinfection of HIV and cholera model using GHF derivative and the solution's behaviour is deciphered. In present work we have used ABC derivative for our proposed model. However, we can extend our work using the GHF derivative for fractional order mathematical models in future.
Despite all the measures taken to control the virus, COVID-19 remained a threat because of its high propagation rate leaving the majority of the population susceptible to the disease and escalating the need to develop a suitable and effective vaccine. The process to gain immunity against COVID-19 is still under discussion. Most of the people infected with COVID-19 developed the immunity against the virus, but still it is unknown how long this immunity will last as there are some reported cases of getting the infection for a second time. In this situation, the only solution to gain herd immunity is mass vaccination. Hence, to achieve this milestone, vaccination is required for a huge portion of the global population. In this regard, WHO helped governments to launch and refine their vaccination programmes along with UNICEF, Gavi, and partners. According to reference [19], 43.5% of the world population has been given at least one shot of COVID-19 vaccine till the 21 st of September 2021. Modified mathematical models with an additional compartment of the vaccinated class are discussed in recent literature to foresee the efficacy, effect, and dynamics of the vaccination. Meng et al. [20] established an SEIRV model to differentiate between the effects of mandatory and voluntary vaccination methods on heterogeneous networks. Ramos et al. [21] investigated the impact of variants and vaccines on the spread of COVID-19 applying the model to Italy as a test case. Tetteh et al. [22] analyzed the mass and ring vaccination strategies to control the transmission of SARS-CoV-2 by varying vaccine efficacy rate. Harizi et al. [23] used the SIRV model to check the future dynamics of the pandemic by varying vaccination rates in Canada. Karabay et al. [24] presented an SEIR simulator to model effective immunization and sterilizing outcomes in order to assess the vaccination effect. Tuteja [25] numerically solved the modified SEIR model including the vaccinated class to investigate the viral dynamics of the epidemic. Ghostine [26] proposed an SEIR model with an extension of the vaccination compartment to investigate the effect of vaccination in Saudi Arabia. Wintachaia and Prathomb [27] presented a detailed stability analysis of the COVID-19 epidemic model to investigate the vaccine efficiency. Alsakaji et al. [28] examined the dynamics of COVID-19 in the UAE using an exptended SEIR epidemic model with vaccination, temporal delays, and random noise. For transmission of the infectious disease, Kou et al. [29] established a multiscale agent-based model to investigate the impact of nonpharmaceutical interventions and vaccination. Acuña-Zegarra. [30] formulated a COVID-19 optimal control problem to study the efficacy and responses of the vaccine induced immunity. Albani et al. [31] estimated the parameters of SEIR like epidemic model using daily reports of COVID-19 of New York and Chicago to analyze the impact caused by vaccination delay.
Motivated from the above studies, we have presented in this work a modified SEIRV epidemic model with an additional compartment of vaccinated individuals to examine the effectiveness and efficacy of the vaccine by varying its rate considering the ABC derivative. We have used ABC fractional order derivative because a number of experimental facts exhibit that the natural dynamics go along with fractional calculus as the fractional order derivatives are nonlocal in behaviour and possess the hereditary properties and memory effects. Also, we have constructed a new, more efficient, and effective numerical scheme for fractional differential equation with Atangana-Baleanu derivative in Caputo sense using Simpson's 1/3 method for the approximation of the integral function involved in the problem.
The paper comprises different sections. In section 2, some basic definitions are recalled which have been used in the work later. SEIRV epidemic model for COVID-19 is presented in section 3. The equilibrium points and the basic reproduction number using next generation method are calculated, and the stability of these points is also analyzed in this section. For the solution of the proposed model modified Simpson's 1/3 method alongside complete stability and convergence analysis, it is presented in section 4. Conversion of mathematical model to fractional order model is given in section 5. Numerical simulations are also presented in the same section.

Preliminaries
This part will review some of the fundamental definitions of fractional calculus, which will be useful in the subsequent work. 2 Computational and Mathematical Methods in Medicine Definition 1. For x ∈ ℍ 1 ð0, TÞ and 0 < ϑ < 1, the ϑth-order ABC derivative and integral of xðtÞ are defined [32] as where JðϑÞ represents the normalization function given by with JðϑÞj ϑ=0,1 = 1 and E ϑ is the Mittag-Leffler function defined as E ϑ ðzÞ = ∑ ∞ n=0 ðz n /Γðnϑ + 1ÞÞ.
The lemma below presents the Gronwall inequality [33], which will be used for the proof of our main result.

Lemma 2.1. The Gronwall inequality for fractional order differential equation states
Let C 1 > 0 independent of δt, C 2 ≥ 0 and {u n } satisfy the inequality with 0 < ϑ ≤ 1: Then where the parameter E ϑ is the Mittage-Leffler function. Particularly when ϑ = 1, inequality (5) becomes In this study, we used an SEIR model with vaccination effectiveness to forecast the COVID-19 condition when a vaccine is released.

Mathematical Model
Mathematical models are of great importance in predicting the behaviour of the viral epidemic diseases. There are a number of models which can be used to describe and investigate the spread of the disease. Recently COVID-19 emerged as a viral epidemic, and we are not aware enough about the behaviour and spread of this virus. To investigate the COVID-19 situation under the effect of vaccination, we consider the following SEIRV epidemic disease model: To define the model equations, the total population, represented by N, is separated into susceptible, exposed, infected, recovered, and vaccinated classes denoted by SðtÞ, EðtÞ, IðtÞ, RðtÞ, and VðtÞ. The flow chart of the above system is shown in Figure 1.
The above system of equation (7) is subject to the initial conditions and the following are the parameters: μ : per-capita natural death rate, γ : the rate at which infectious people recover, ρ : average fatality rate for virus-infected people, β : disease transmission probability per contact times the amount of contacts made per unit of time (dimensionless), Λ : the birth rate per-capita, ε : progression rate from E to I, 1/γ : contagious period, 1/ε : incubation period, δ : vaccination rate, σ : vaccine ineffectiveness, 0 ≤ σ ≤ 1, 1 − σ : vaccine effectiveness, having T units per unit of time.

Stability Analysis of SEIRV COVID-19
Model. In this section, stability analysis will be carried out to find the disease-free equilibrium and endemic equilibrium points. In order to get these points, each of the equations given in system of equation (7) will be equated to zero.

Disease-Free Equilibrium and the Basic Reproduction
Number ðR 0 Þ. Because there is no disease propagation in disease-free equilibrium, so I = 0, and hence we get so the disease-free equilibrium point is E 0 ðS 0 , E 0 , I 0 , R 0 , V 0 Þ = ðΛ/ðμ + δÞ, 0, 0, 0, δΛ/μðμ + δÞÞ. The average number of secondary infections induced by an infectious individual is given by the basic reproduction number R 0 . It estimates the growth of the virus outbreak. The basic reproduction number R 0 is determined using the 3 Computational and Mathematical Methods in Medicine next generation method [34]. Using _ E and _ I from system of equation (7), the Jacobian matrices of the disease-free equilibrium are as follows: Here FU −1 is the next generation matrix and we can obtain the basic reproduction number as R 0 = ρðFU −1 Þ, so Substituting the values of S 0 and V 0 , the following value of the reproduction number is obtained for the system of equation (7): Theorem 3.1. If R 0 < 1 then the COVID-19 model of equation (7) disease-free equilibrium point is locally asymptotically stable.

Existence and Uniqueness of Endemic Equilibrium.
In order to investigate the existence and uniqueness of the endemic equilibrium E * for R 0 > 1, we will set all the derivatives of model of equation (7) equals to zero. Let E * = ðS * , E * , I * , R * , V * Þ be the endemic equilibrium of model of equation (7). Then by solving the following equations we can obtain the endemic equilibrium.
Now, solving first, third, and fourth equations of equation (17), we get 4 Computational and Mathematical Methods in Medicine Also equation (17) gives The following quadratic equation for I * is obtained: where It is evident that A is always positive while C is negative for R 0 > 1. Hence, there exists a positive and unique value of I * , which results in a unique endemic equilibrium for R 0 > 1.

Numerical Approximation
Consider the nonlinear fractional order differential equation where 0 < ϑ < 1: The vector function f in the differential equation is real valued and continuous and satisfies the Lipschitz condition.
To derive the numerical scheme for the above mentioned nonlinear fractional order differential equation, the given interval ½0, T will be subdivided into equal intervals with time step δτ = T/N, where δτ n+1 − δτ n = δτ for n = 0, 1, 2, ⋯, N − 1: For convenience, the approximate solution of xðτ n Þ will be denoted by x n with δτ n = nδτ, n = 0, 1, 2, ⋯, N.
Applying the integral operator of equation (2) on equation (22), we get where the integral I ϑ in the above equation is the ABC inte-gral. Using equation (2) to discretize the integral I ϑ , we get f ðxðsÞÞ is approximated on each subinterval ½τ i , τ i+1 , by the following piecewise quadratic interpolation polynomial.

Lemma 4.2. For
Proof. According to the definition of b i ði = 1, 2, ⋯, n − 2Þ, it deduces that Computational and Mathematical Methods in Medicine implies Note that 0 < ϑ < 1, so the coefficients of power of 1/ðn − iÞ are less than 1. Thus, the series given above is convergent, and hence the proof is completed for b i ði = 1, 2, ⋯, n − 2Þ. Similar strategy can be adopted for b n , so we neglect the proof.

Numerical Results and Discussion
Considering the case of Italy, we applied the proposed model to analyse the impact of vaccination. The vaccination campaign for COVID-19 started on 27 th December 2020 in Italy managed by the Ministry of Health. In the early months of the vaccination campaign, the Italian government targeted health doctors and administrative staff as well as guests and nursing home personnel followed the vaccination of elderly people and public service personnel in second phase of the campaign. Italian government bought four out of seven, WHO approved vaccines to facilitate its people against the virus. So far, 73.3% of the total population has received at least one dose and 65.5% of the total population is fully vaccinated. As of 17 th September 2021, top ten

10
Computational and Mathematical Methods in Medicine regions of Italy for the total number of administered COVID-19 vaccine doses can be seen in Figure 2. Moreover the number of people who received at least one dose and fully vaccinated according to the data [19] are presented in Figure 3. Further, graphical presentation of the numerical simulation for the fractional model (51) are given in this section. Table 1 shows the estimated model parameters and initial values used to evaluate the impact of vaccination on the COVID-19 epidemic. Figure 4 shows the impact by varying the vaccination rates on different classes of the proposed model with δ = 0 representing the model output without vaccination. Although there is no remarkable change in the graph of susceptible class, we can observe a significant decrease in infected and exposed classes by increasing the vaccination rate. We can see the peak of the infected as well as that of exposed individuals gradually decreases with the increase in effective rate of vaccination. From Figure 4, we can observe that for exposed class there is a difference of almost 7 × 10 6 individuals in peak values for no vaccination to a vaccination rate δ = 0:012. Similarly, the peak values for no vaccination compared to a vaccination rate of δ = 0:012 for the infected class differ by almost 7 × 10 6 individuals. So, interpreting graphically the dynamics of the proposed model to observe the effect of vaccination by altering the vaccination rate suggest that, mass vaccination along with the implementation of the control measures until the achievement of complete herd immunity is needed to control the spread of the disease. Hence, majority of susceptible population should be vaccinated to prevent the new waves of COVID-19.  Figure 4: Impact of vaccination rate on all the four susceptible, exposed, infected, and recovered classes.

Computational and Mathematical Methods in Medicine
The vaccination of the susceptible class will control the spread of the disease resulting in reduction of the mortality rate, which can clearly be seen by recovered class presented in Figure 4. Simulation results show that immunization is an essential strategy for lowering the incidence and severity of COVID-19 infections in the workplace, in communities, and across the globe. Therefore, it is advocated for policymakers to increase public trust and confidence in the efficiency and security of vaccinations as well as the expertise and dependability of the institutions that provide them. Figures 5 and 6 show a comparison of integer order model, fractional model for fractional order ϑ = 0:975,0:95,0:925; ;0:9, and real data of total number of deaths D = ϑI owing to the COVID-19 outbreak in Italy. Figure 5 presents the comparison for first 80 days starting from 23 rd February 2020 to 13 th May 2020 of the pandemic before vaccination. We can see that after 40 days, real data is better explained by fractional order model in comparison to integer order. Figure 6 is assigned to display the comparison from 1 st July 2021 to 19 th September 2021 along with administered vaccination effect considering the initial values of the compartments as of 1 st July 2021 . The relevant COVID-19 data can be found at https://github.com/pcmdpc/COVID-19. Also we can see from our proposed model, fractional order models are better at fitting the real data than classical models because fractional order derivative and integral operators, due to their nonlocal nature, depend not only on their current state but also on all of their previous states, making them more effective than other classical deterministic operators in predicting the biological model's future state.

Concluding Remarks
Within the Atangana-Baleanu Caputo fractional operator, we studied alternative conditions for qualitative characterization and dynamics of the suggested model. Calculating disease-free and endemic equilibrium points, as well as the reproduction number, are used to examine the prerequisites for local stability. A newly developed numerical scheme is used to obtain the solution of the concerned fractional order model. The proposed scheme's stability and convergence are briefly examined. To depict the behaviour of the COVID-19 virus, a graphical description of the dynamics of the relevant system is developed within different fractional orders under the effect of vaccination suggesting the need of mass vaccination to gain the herd immunity in order to prevent the further spread of the disease. Due to the limited access of the available vaccine, other preventive measures should also be continued.

Data Availability
The data used to support the findings of this study can be found at https://ourworldindata.org/covid-vaccinations.

Conflicts of Interest
The authors declare that they have no conflicts of interest.

1.305
Dead individuals ×10 5 Figure 6: Comparison of the number of people who died for various fractional orders vs. real data (after vaccination).