Dynamics of SEAIQR Model with Saturated type Treatment: A case Study of Spain COVID-19

Background. Outbreak of the Covid-19 is now an ongoing global health emergency. At the end of December 2019, the first infection was reported in Wuhan and the world did not pay attention to this extremely contaminated disease and plucked to react rapidly. The World is in an vulnerable state in disease spreading, facing a great loss of lives and socioeconomic aspects also. That is why we have proposed a potential mathematical model with data analysis to predict and control the outcome of this pandemic. Methods. The model presented the epidemic dynamics of multiple compartments. We collected available online data of Spain. In primary step, we estimated the parameters using either the data analysis or reference papers. Then we did the data fitting analysis in comparison with the outcome of our mathematical results. The results of the system depended not only the parameters also on social consciousness. Results. It is found that disease progression in this model is determined by the basic reproductive ratio, R0, the actual epidemic of R0 and effective R(t) of each day. If R0 > 1, the number of latently infected individuals grows exponentially; endemic solution is stable while infection rate decays if R0 < 1. The optimal control theory stated that vaccination and treatment strategies are highly effective for reducing both susceptible and infected population and to increase the recover rate high. In Spain, after state of alarm (quarantine) on 14 March 2020, reported cases increasing for 13 days only and from the 14th day, daily reported cases started to decline albeit with small fluctuation. Our proposed model approximates that the disease in Spain could be fully under control by after July 2020. Conclusion. Outbreak will be in control of health care system, reduce the death rate and will ensure social-economic stability.


Introduction
In COVID-19 pandemic, the whole world were apparently unprepared or might did not pay well attention at the earlier stage after spreading in Wuhan, China. The disease caused by SARS-CoV-2 is now having a devastating impact across the world. At the current stage, almost all countries and territories are a victim with 4,269,782 reported cases with 287,529 deaths on May 12, 2020 and this is a global health emergency [1,2]. Compared to reported cases, the recovery rate is still very slow.
The outbreak of severe acute respiratory syndrome SARS, appeared in November 2002 in China. A newly emerged coronavirus was responsible for this outbreak which is known as SARS-CoV [3]. Various experimental data indicated that wild animals were the source of the SARS. Serological experiment of the people associated with those wild animals revealed that the animal traders consisted larger amount of IgG, an antibody active against the SARS-CoV compared to those who were the traders of vegetables [4]. Evidence showed that 3 out of 4 patients were supposed to have either straight or secondary interaction with the masked palm civets. Among these three patients, one of them was waitress in TDLR restaurant where the animal was kept for bizarre dishes, another one took his meal in that restaurant and the dining table was close to the cage where the animal was kept and the third patient took his meal in another adjacent restaurant. It is remarked that this study did not say anything about the natural reservoir of the virus [5,6]. Later study revealed that several species of horseshoe bats were the natural reservoir of SARS like-CoV. Among four species, R. pearsoni and R. macrotis showed affirmative data in serological and PCR test and the remaining species R. pussilus and R. ferrumequinum showed a positive result in either PCR or serological test [7].
Middle East respiratory syndrome coronavirus MERS-CoV took place in Saudi Arabia in 2012 where dromedary camels were thought to be the intermediate source for the transmission of the virus [8]. Research data suggested that about 55% of primary victims of MERS-CoV had direct contact with the dromedary camels or foods of camels and the remaining 45% portion of the victims did not have any evidence of contact with the camels or infected individuals. Hence still its a question, how they got infected by the virus was unknown [9]. Camels were declared as the zoonotic host for human infection, evidence suggested that bats were the natural reservoir host of the virus [10]. Small genomic sequence analysis of bat acknowledged bat as the natural reservoir of MERS-CoV [11].
As a continuation, the virus SARS-CoV-2 was related to the coronavirus responsible for the SARS outbreak of 2003, and the transmission of the virus was also zoonotic [12]. In March 13, 2020, the World Health Organization considered Europe as hot spot of the 2019-20 coronavirus pandemic, and by April 24 (2020), the maximum part of the Glove were affected by the outbreak [1, 2,13]. The daily reported cases were being doubled over periods of 2 to 4 days by country across Europe, North-America, Asia [14], and specially the phenomenon was drawn a lot of attention for these zones.
The present hypothesis about SARS-CoV-2 is that the outbreak started in bats, then moved to another unidentified species before being transmitted to humans [6,15]. Additionally, in recent time COVID-19 reinfection is drawing people's attention. The World Health Organization (WHO) investigating reports of COVID-19 reinfection and the organization is working with clinical experts to get more information on this matter [16]. In South Korea 91 patients who had tested positive, and declared negative after two consecutive tests, tested positive again while 14% of recovered coronavirus patients in China have tested positive again [16,17].
Since the world is dealing with a highly contagious disease without medication, right now the only effectual way for us to protect is prevention. Social and physical distancing, lock-down and testing are some measures prescribed by WHO to control the outbreak. Data analysis and mathematical modelling are one of the core component to study these undertaken polices for the best outcome. The mathematical growth map of COVID-19 was introduced using two well known growth functions and establised that the growth is exponential [18]. In mathematical modelling, some recent studies provided different guidelines by introducing basic reproduction number, education and socio-economic index and lockdown strategies [19,20,21,22,23]. A recent study developed a COVID-19 epidemic model with five compartments which analyzed both theoretically and numerically by evaluating the equilibrium points and their stability analysis. The dynamics was shown for a different group of the population over time and the findings was based on the sensitivity analysis; highlighted the role of outbreak of the virus that might be useful to avoid a massive collapse a developing country, Bangladesh [22]. In the literature, an ordinary differential equation (ODE) model studied to drag down the spread of COVID-19, with the influnce of social consciousness towards the epidemic; where the prescribed model determines the ongoing social response project the scenario of COVID-19 handled with a higher social consciousness level [23]. The dynamics of SEIR model and its optimal policy to control the disease in Italy were studied in [24].
For any data-driven epidemic models, the basic reproduction number, R 0 , is an important measure of the strength of any pandemics. In susceptible population, R 0 is defined as the number of secondary infections generated by one primary case. For any outbreak of disease, if the estimate of R 0 is in early stage then the guidance of control measures can be significant. To calculate the basic reproductive number, we have two phased (i) R 0 for the actual epidemics estimated by parameters, and (ii) the analysis of it's time evolution, R(t), since it changes continuously during the outbreaks [25,26].
In this paper, we have studied the current scenario of the world and the European country Spain, as a special case since this territory have higher socio-economic index [13,27,28]. The study presents an analysis of various control mechanisms using mathematical model and abridged data fitting based on ongoing viral epidemic of the world. The key objectives of this study are as follows: We will work with real-life available discrete data of Spain to understand the cases and project the control of the infection.
We expect that the prediction of controlling measure of COVID-19 will be able to validate the dynamics of SEAQIR model to get better accuracy of results.
The main novelties in this study are 1. Analytical and Theoretical results are established and presented in terms of reproductive ratio, R 0 . We have obtained the stable DFE for R 0 < 1 and endemic solution as long as R 0 > 1.
2. Using central manifold theory, we have established the results for backward bifurcation while R 0 ≡ 1 with some other conditions.
3. The basic reproduction number for actual epidemic, R 0 and time evolution outbreak, R(t) are calculated and their respective results are established.
4. Numerical demonstrations validated the analytic outcomes.
5. Using optimal control theory, it is concludes that vaccination and treatment are highly effective and gives more number of recovered population. 6. Data analysis and model results predict the time boundary to control the epidemic in Spain.
7. Besides data analysis, the proposed quarantine model shows that human population can control and protect the spread of infectious diseases by creating social isolation, hospitalization, lock-down and physical distancing.
The paper is organized as follows: Mathematical Model is elaborately discussed in Section 2. Positivity and boundedness of solution including auxiliary results are described in subsection 3. Stability analysis and bifurcation analysis are prescribed in Section 4. Section 5 is accomplished the data analysis in comparison with model solution using parameter estimation and model validation with further prediction to control the endemic, as a case study in Spain. We have studied the optimal control theory for vaccination and treatment strategies to increase the recovered rate in Section 6. Finally, Section 7 outlines the summary and discussion of the results.
In the following section, we will discuss our mathematical model and formulation of the model elaborately.

Model Formulation
Scientifically it is proven and well established that mathematical and computational models can help to understand biological scenarios and can predict epidemiological aspects. In Infection rate of E class σ 1 Recovery rate of I class the proposed model, let the total population is divided into six compartment those are the susceptible population (S), exposed population (E), asymptotic population (A), quarantine population (Q), infected population (I), and recovered population (R). We assume that the disease COVID-19 spreads mainly from A and I compartments and the rate of incidence is of the from (αA + βI)S. Since a large number of persons is infecting by this disease in very short time but the number of hospital is limited. Due to limited medical facility we consider the rate of treatment is non-linear and of the from aI 1+bI . The flow diagram of the proposed model is given in Figure 1 and the corresponding mathematical equations are given in (2.1).
where the interpretation of parameters is presented in Table 1.
In the next section we shall establish the positivity and boundedness of solution of (2.1). Positivity of solution in biological model is most important because it establish the non-negativity of the solution curve i.e. non-negativity of the considered compartments.

Positivity and Boundedness of Solutions
Consider the first equation of (2.1) and multiply it by integrating factor exp(αA(t)+βI(t)+ ), integrating and using the initial condition we get It is obvious from the above expression S(t) ≥ 0. Using the similar arguments it can be easily shown that other variables are also non-negative.
To prove the boundedness of the solutions add the equations in (2.1) and using the notation for total population, N (t) = S(t) + E(t) + A(t) + Q(t) + I(t) + R(t), we have Integrating and using the initial conditions and lim sup for t → ∞, we get Thus we can summarize the above results in the following theorem: Theorem 1. All the solutions are feasible and the set Ω = (S, E, A, Q, I, R) ∈ R 6 + : S + E + A + Q + I + R ≤ Λ is a positively invariant set for the system (2.1).
In the next section we shall establish the existence of equilibrium points, the expressions of the basic reproduction number (R 0 ).

Stability and Bifurcation Analysis
In this section we first find the disease free equilibrium point (DFE) and then using the next generation matrix approach the basic reproduction number will be found. The stability of the equilibrium point will be discussed in terms the basic reproductive ratio. Then backward bifurcation will be discussed.

Basic reproduction number, DFE and EEP
The basic reproduction number is denoted by R 0 is defined as expected number of secondary infections that a single infected individual generates [29,30]. In epidemic system it has great importance as it is the indicator of persistence or irradiation of disease. Using next generation matrix [30,31] the basic reproduction of (2.1) is found here. Since the DFE is E 0 λ , 0, 0, 0, 0, 0 and hence the basic reproduction number can be found using the analytical approach. Let us we define F = ∂F ∂x j (E 0 ) and V = ∂V ∂x j (E 0 ) then the next generation matrix F V −1 and the spectral radius of (F V −1 )=R 0 is the basic reproduction number is given below Let E * (S * , E * , A * , Q * , I * , R * ) be the endemic equilibrium point (EEP) of the considered system then it satisfies the steady state condition for the system (2.1). Solving the steady equations we obtain 1+bI * + σ 1 I * and E * satisfies the equation It is clear from the expression of coefficients the equation (4.1) may have either two positive roots or no positive root if R 0 < 1. For R 0 > 1 the number of positive roots is one or three if the positivity condition of J * is satisfied. Thus when R 0 < 1 the model (2.1) have at most two endemic equilibrium points and for R 0 > 1 the number of equilibrium points are three or one.

Stability of disease-free equilibrium state (E 0 )
To eradication of disease from the system the stability of the DFE is essential.
Theorem 2. The DFE will be locally asymptotically stable if R 0 < 1 and unstable if R 0 > 1.
Case I: Suppose R 0 > 1. Then P (0) > 0. Again, lim λ→∞ P (λ) = −1. Here P (λ) is a continuous function of λ and so by applying Bolzano theorem on continuous function we have P (λ k ) = 0 for some λ k > 0. Thus, at least one eigenvalue of the Jacobian matrix J(E 0 ) must be positive and so the DFE, E 0 is unstable.
Case II: Suppose R 0 < 1 which yields P (0) < 0. If possible let us assume that P (λ) = 0 has at least one root of the form x + iy, where x ≥ 0 and x, y ∈ R. Then P (x + iy) = 0.
Hence, all the roots of the equation P (λ) = 0 have the form x + iy, where x < 0 and x, y ∈ R. Thus, in this case the DFE is locally asymptotically stable. Hence the theorem is proved.

Backward Bifurcation
In this section, we shall investigate the existence of backward bifurcation of the considered COVID-19 model. Biologically, this type of bifurcation is most important because usually the disease eradication happens when R 0 < 1 but in this situation there exists another stable endemic equilibrium point for R 0 < 1 and so eradication of the disease depends not only on the value of basic reproduction number but also on the initial number of infection level. When R 0 = 1 i.e. when then one of the eigenvalues of the Jacobian matrix corresponding to the system (2.1) at DFE is zero. Using the Theorem by Castillo-Chavez and Song [29,31,33], we shall establish the existence of backward bifurcation of the system (2.1).
Theorem 3. Suppose c 3 > 0 and c 1 c 2 > c 3 . Then the system (2.1) undergoes the backward bifurcation at R 0 = 1 with respect to the parameter α if φ > 0, where c 1 , c 2 , c 3 and φ are defined in the text.
Proof. First we set At the DFE E 0 the Jacobian matrix J(E 0 )| α=α [BB] has the eigenvalues 0, − , − and other three are the roots of the following cubic: It is obvious that c 1 > 0 and R 0 = 1 implies c 2 > 0. Hence by the given conditions we can conclude that other three eigenvalues have negative real parts. Thus, we can apply the Theorem by Castillo-Chavez and Song to determine the direction of the bifurcation. Let W and V be the right and left eigenvector of J(E 0 )| α=α [BB] corresponding to the zero eigenvalue then and the coefficient Since the coefficient ψ is obviously a positive number, hence the system (2.1) exhibits the backward bifurcation if φ > 0. Hence the theorem is proved.
It is clear from Figure 2 that if we vary the transmission rate of infection (α) from the asymptomatic class, then the system generates two endemic equilibrium points for R * 0 < R 0 < 1, where R * 0 is some positive quantity, among them one with lower infection level is unstable along with the stable DFE and other with higher infection level is stable. Epidemically, the basic reproduction number less than unity is not enough to eradicate the disease completely if the conditions for backward bifurcation hold.

Stability of endemic equilibrium state
The Jacobian corresponding to the endemic equilibrium point E * is J(E * ) is given below One characteristic root of the above matrix is − and other five satisfies the following equation Since one root of the characteristic equation is negative and the system will be locally asymptotically stable if other roots are negative or real part is negative. The other roots will be negative or negative real part if the Routh-Hurwitz [29] is satisfied i.e. a 1 , a 2 , a 3 , a 4 , a 5 > 0, a 1 a 2 a 3 > a 3 2 + a 1 2 a 4 and (a 1 a 4 − a 5 )(a 1 a 2 a 3 − a 3 2 − a 1 2 a 4 ) > a 5 (a 1 a 2 − a 3 ) 2 + a 1 a 5 2 . Since expressions of the coefficients a i s, i = 1 − 5 are large, verifying the above conditions is difficult task. Numerically we shall justify the condition of stability of the endemic equilibrium point.
Since Z(t i , θ) is function of the model parameters θ ∈ Θ, our aim in the present section is to find θ ∈ Θ such that is minimized subject to θ ∈ Θ, n is the number of observations [34]. There are several optimization technique to minimize the error for suitable θ ∈ Θ, here we have used the Matlab package fminsearch to estimate important model parameters. Rest other parameters are collected from published research documents.

Sensitivity index
The sensitivity analysis reveals parametric influence of the disease spreading in epidemic modelling. Usually the used methods for sensitivity analysis are the perturbation of fixed point analysis and the uncertainty in parameter estimation. Using the sensitivity index one can estimate the rate of change of variables when the parameter changes [35,36]. On the other hand the most important tool to study the dynamics of epidemic model is R 0 and it is function of the model parameters. Here our goal is to determine the significant model parameters which are controlling R 0 . To study the effect of sensitivity here we use global sensitivity index [35] for R 0 with respect to the model parameter φ is denoted by Γ φ R 0 and is defined by If Γ φ R 0 = ξ then y% increase (decrease) of φ will increase (decrease) ξy% when ξ is positive and for negative ξ increase of parameter will decrease R 0 .

Confidence interval of the model Parameters
The confidence interval for the model parameter can be computed using the sensitivity matrix as described in [34] and used by [37]. The sensitivity matrix (J) for the parameters Θ(θ 1 , θ 2 , ..., θ l ) is defined as  Using the complex-step approximation method the partial derivative ∂Zt k ∂θ j can be calculated [38,39]. Using complex step ih for small h (for computation purpose we use h = 10 −40 , i is the unit complex number). Using Tayor's The imaginary part of the both side gives the corresponding partial derivative where O(h 2 ) represents the second and higher order terms. From the above J, one can be obtained then the standard order (SE θ j * ) for the parameter θ j as defined is computed for the basic estimated parameter (θ * ). Using the above results the confidence interval can be obtained.
Since we have data availability only for cumulative infected populations (Z-class). Using the method defined above we found that confidence interval for most of the model parameters are negative and the SSE is large. The uncertainty embedded in the large interval can be minimize fixing in the less effective model parameters. This processes was successfully applied in the previous studies [37] and not violated the significance of data fitting. For this purpose we first compute the sensitivity matrix and the most sensitive parameters are presented in the box diagram (see Figure 3). We have computed the confidence interval for the most sensitive model parameters and other parameters are fixed as initial least square estimation level.

Data management
In this study we considered the country Spain drastically infected by SARS-CoV-2 virus and collected available online data in chronological order from Worldometer [1, 2]. We considered Spain since it was one of the most affected country in the world at the time of this study. We included the data from February 14 to April 28, 2020 (total 65 days).
To ensure the curve fitting we did not include any additional data outside the time interval as declared above. For simplicity and to compare the results, zero (0) was used to denote the first day when a case was reported for Spain; see graphical presentations in later section.

Case study: Spain
To study the particular case: the spreading of COVID-19 in Spain we have considered the real cases in cumulative number of infection (from Feb 14 to April 28 2020) and using the above described formula we have estimated the model parameters, using the global sensitivity index method we found the sensitivity index of the parameters. To calibrate the numerical simulation the initial population size of the susceptible and infected population is taken as S(0) = 46751619 and I(0) = 3 and other four initial conditions are estimated as E(0) = 10, Q(0) = 1, A(0) = 2, and R(0) = 0. The estimated model parameters, sensitivity index and the corresponding sources are in Table 2. Using the procedure described in section 5.3 we have identified the most sensitive model parameters which is given in the box plot of the parameters and the confidence interval for the parameters, see Figure 3.
It is clear from the Table 2 that the first six parameters have positive correlation with R 0 and the rest nine parameters have negative correlation with R 0 . With the increase of positive correlated model parameters the prevalence of infection will increase and for negatively related parameters will decrease the number of infection. The best model predicted cumulative curve and the real cumulative, bar presentation of the daily infected cases & the model predicted values and corresponding residuals is given in Figure 4. It is clear from Figure 4 (c) that the residuals are randomly distributed. The random distribution of the residuals imply that the fitness of the data with the model is over all good [29].

Model Validation
In previous section we consider the real data of Spain to predict the model parameters up to April 28, 2020. To validate the model we have presented the model predicted per day new infection and the real cases in Figure 5 from April 28 to May 10, 2020. It is clear from the figure there are minimum variation of the model prediction and the real cases.

Effect of Significant model parameters
Now we shall investigate the effect of significant model parameters in the spreading of COVID-19 in Spain. Keeping other parameters fixed as in second column of Table 2 we varied the parameters α, β, µ 1 γ 1 and γ 4 once at a time.
(a) Effect of α, contact rate between S and A class In Figure 6 we present the time series of the infected population for different values of the model parameter α in the estimated interval as given in Table 2. It is clear from the sensitivity index that with the increase of α the basic reproduction number increases i.e. prevalence of the disease will increase. In the figure the red curve denotes time series of the infected population for the model predicted parametric values, the upper and lower curve denotes the same for upper limit and lower limit of α in the interval, the space between the curves denotes the prevalence of infection for different values of infection rate between the asymptotic and susceptible class in the estimated interval. It is clear from the figure that for higher values of α the daily new infection will reach the highest value shortly then it will start to decrease.   To minimize value of µ 1 we have to quarantine the exposed properly then spreading of disease can be controlled easily.
(d) Effect of γ 1 , transmission rate of E class to A Since γ 1 is the rate at which the exposed class transfer to the asymptotic class. If γ 1 higher then a large number of infectious people will enter the asymptotic compartment and consequently the disease will spread among the susceptible people shortly. For lower values of γ 1 disease spreading will be lower. In Figure 9 we have presented the time series of the infected persons for different values of γ 1 . The description of the graph is similar as the last three. It is clear from the figure to control the disease spreading we have to control the parameter γ 1 . (e) Effect of γ 4 , recovery rate of E class In Figure 10 we have presented the time series of the infected population for different values of recovery rate of the exposed class (γ 4 ) as given in the confidence interval. Like the previous two figures here the red line is the model predicted daily new infection for COVID-19 but the lower line is for highest value of γ 4 and the upper line is for lowest value. Hence with the increase of this parameter the number of infected population decreases i.e. to control the disease the recovery rate of the infected population plays important role. The recover rate will increase if the immunity of the infected person is high. So, to get recovery form the disease individual have to gain high immunity through life style.

Contribution of R 0 by infectious and asymptotic persons
It is clear from the expression of R 0 , the two parts are respectively arises due to the effect of asymptotic and infectious class. Using the estimated parameter values we obtain R 0 = 6.692558249 and the value of two parts are R 01 = 5.296117140 and R 02 = 1.396441109. It is clear from the values that the asymptotic class gives the 79% contribution and the infectious class gives 21% contribution in the basic reproduction number. Thus the asymptotic class gives most contribution i.e. in spreading the COVID-19 infection, the asymptotic class take most important part compare to the infectious individuals. Thus to control the disease of we have to minimize the asymptotic persons and after test we have to transfer them in the infectious compartment.

Estimation of R 0 for actual epidemic
Several methods are existing to estimate the reproductive ratio, R 0 for vector-transmitted diseases using actual data of epidemics. In this paper, we calculated R 0 from the initial growth phase of the diseases. At the primary stage of the epidemic, we assumed that the number of cumulative cases, H(t), varies as H = exp(ζt), where ζ is the force of infection. In this case, the number of asymptotic, exposed and infectious vectors vary similarly, where A 0 , E 0 , Q 0 , and I 0 are constants. Also form the first equation of (2.1), the constant susceptible population is S 0 = Λ . Assume that delay in treatment parameter, b = 0. Moreover, the number of recovered vectors (non-susceptible class) can be assumed negligible such that Substituting equations (5.2) and (5.3) in the respective expression for the derivative as presented in the system (2.1), we obtain first form E−equation such that Similarly, we obtain from A, Q and I− equations, respectively For simplicity, using the notations D 0 , D 1 , D 2 , and D 3 as defined in section 4.1, we have Using the above four relation the association of R 0 and the force of infection can be done in the form where P 0 = ζ D 0 + 1, P 1 = ζ D 1 + 1, P 2 = ζ D 2 + 1 and P 3 = ζ D 3 +a + 1. Now it is time to estimate the force of infection, ζ. Since the number of new cases per day, I(t), is the derivative of I(t) in relation to time t, then at the beginning of the epidemics, I(t) ∼ I 0 exp(ζt). Plotting the number of new cases per day against the cumulative number of cases Z(t), the phase of exponential growth of the cumulative number of cases is evidenced by a linear growth of the curve, the slope of which is the force of infection and which can be computed by a least-square linear fit of this linear phase [41,42]. For the data of Spain in 14 Feb to 28 April, 2020 for the COVID-19 outbreak, shown in Figure 11(a), we obtain the calculated range of ζ = 0.1497553492834 ± 0.0091 day −1 based on the slope shown in Figure 11(b). The straight line (red color) indicates the growing linear parts of the plots corresponding to the initial exponential growth of the epidemics.
Substituting the parameter values as presented in Table 2 (no systematic vector control was in place during this epidemic) in expression (5.5), we obtain an estimation of R 0 = 5.4925 and the range is R 0 ∈ (5.1514, 5.8448).

Effective reproduction number
In epidemiology the basic reproduction number (R 0 ) denotes the average number of secondary infection produced by one infected individual during the course of mean infection time. But usually when disease start to spread initially the secondary number of infected in higher and it reaches to a peak point then it decrease i.e. the reproduction number is not always fixed. Our aim in this part to define the time varying reproduction i.e. the reproduction number per day for COVID-19 and determine it using the available method. This type of estimation of the reproduction number is defined as the effective reproduction number at time t and denoted by R(t) [41,42,43]. Here R(t) yields the information about the invade of infection among the susceptible populations with time. To adopt the control policy of COVID-19 in Spain for day to day it is most important to know the time varying reproduction number. The authors in [41,42] defined the renewal equation in the form where b(t) is the number of new cases at t th day and g(τ ) is the generation interval distribution for the disease. Using the procedure described in [41,42] to the model (2.1) the generation distribution can be obtained. Let the rate of leaving of the exposed, asymptotic people are respectively p 1 = + γ 1 + γ 2 + γ 3 + γ 4 , p 2 = δ 1 + δ 2 + and p 3 = d + + σ 1 the generation interval distribution can be obtained as the combination of p 1 e −p 1 t , p 2 e −p 2 t and p 1 e −p 1 t in the form ,j =i (p j − p i ) and mean of the distribution is T = 1 p 1 + 1 p 2 + 1 p 3 and τ > 0. The above relation valid when the force of infection ζ > min {−p 1 , −p 2 , −p 3 }. From the expression of F (t) the effective reproduction per day for the COVID-19 model can be estimated as R(t) = F (t). The value of R(t) for the present model for the date 14 Feb to 28 April, 2020 of Spain COVID-19 spreading is given in Figure 12. The value of R(t) strictly greater than one upto 3 rd April 2020 and then it move up and down of one. Thus growth of disease start to decrease 3 rd April.

Optimal Control
In this section, we have formulated an optimal control problem by introducing time dependent controls due to lock down adapted by some countries to maintain social distancing, causing the reduction of transmission rate of infections. Again, for applying lock down policy by the countries they face huge financial losses and for this reason some countries are adapting partial lock down. On the other hand, in spite of huge financial losses some countries are adapting total lock down to control the spreading of COVID-19. The control of spreading of COVID-19 infections depends on how the lock down is maintained. For the implementation of lock down, we can replace the incidence rate αAS + βIS in the system (2.1) by α(1 − u 1 )AS + β(1 − u 2 )IS, where two controls u 1 and u 2 , respectively, measure the reduction of transmission from asymptomatic and infectious individuals due to social distancing. Thus we formulate an optimal control problem as follows: Here, we construct the cost functional as follows: where the constants B 1 , B 2 , B 3 , B 4 are, respectively, the losses due to presence of exposed, asymptomatic, quarantined and infectious individuals. The constants B 5 and B 6 are the losses associated with respective controls. Main objective of this study is to minimize the COVID-19 infections as well as the financial losses due to implement of lock down for the time T . In order to obtain optimal control strategies we have to find optimal functions u * 1 (t) and u * 2 (t) such that J(u * 1 , u * 2 ) = min{J(u 1 , u 2 ), (u 1 , u 2 ) ∈ U }, where the control set U = {(u 1 , u 2 )/u i (t) is Lebesgue measurable on [0, 1] and 0 ≤ u 1 (t), u 2 (t) ≤ 1 for all t ∈ [0, T ]}. Now, the existence of such optimal control functions for the system (6.1) has been shown by the following theorem.
Theorem 4. There exists optimal control functions u * 1 (t) and u * 2 (t) that minimize J over U .

Conclusion
The time varying pattern mathematical model SEAQIR of the pandemic COVID-19 presented several interesting features in this paper. Considering all major parameters of the progression of the disease, several analytic results established. All the properties necessary for epidemiological relevance have also been proved. We estimated the parametric values for Spain using the existing online available data. The acceptable agreement between data analysis and numerical solutions were established. Theoretical results established and presented in terms of basic reproduction number and the dynamics of the model are determined by the threshold levels of R 0 . We have disease free solution for R 0 < 1 and endemic solution while R 0 > 1. Both theoretically and numerically with data analysis we demonstrated that the infected population will be decrease asymptotically as long as R 0 less than 1. Data analysis and model results predicted the time boundary to control the epidemic in Spain. Besides data analysis, the proposed model showed the vaccination and treatment strategies to control the epidemic; human population can control and protect the spread of infectious diseases by creating social isolation, hospitalization, lock-down and physical distancing. The model ensured the efficiency of the approach to intercept the disease by the limitation of contacts between the individuals through quarantine of infected individuals.
The proposed model confirmed the efficiency of the approach to reduce and gradually stop the disease by the limitation of contacts between the individuals through quarantine of infected individuals. In the model results, we compared outbreak data from Spain with the quarantine model to calculate the growth and decay of the current outbreak. In Spain, lock-down (quarantine) was initiated on 14 March-2020. The maximum number of reported cases was observed on 26 March-2020; which means that the number of infections was on an upward trend for only 13 days. After that day, the number of daily reported cases was observed to decrease asymptotically after which the health authorities appeared to gain control of the situation, which is an the most important interpretation of mathematical model. These values can be used by policy makers and other professionals to make decisions on how certain interventions such as distancing and isolation can be effectively implemented. Similar scenarios can be observed in several European countries such as Germany, UK and France.
The current viral epidemic pushed every medical system to its breaking point without any permanent solution to the disease. So far, there has been minimal advancement towards the use of medication or vaccination to control the spread of COVID-19. As such the most effective mitigation measure currently available is physical distancing to reduce interaction between infected individuals and susceptible populations. The efficiency of quarantine methods determines our protection against the virus.
In this study there were one notable limitation, for example: the model taken into account the beginning of the disease development and we ignored the spatial distribution (space x) of population densities. Instead we considered the case when infection rate and other compartments population were changing with time. Due to these limitations, the accuracy of the existing results and prediction in this paper could vary for a short period