# Global Analysis and Optimal Control Model of COVID-19

Authors: Sacrifice Nana-Kyere; Francis Agyei Boateng; Paddy Jonathan; Anthony Donkor; Glory Kofi Hoggar; Banon Desmond Titus; Daniel Kwarteng; Isaac Kwasi Adu
Journal: Computational and Mathematical Methods in Medicine
 (2022)

---

## Abstract

COVID-19 remains the concern of the globe as governments struggle to defeat the pandemic. Understanding the dynamics of the epidemic is as important as detecting and treatment of infected individuals. Mathematical models play a crucial role in exploring the dynamics of the outbreak by deducing strategies paramount for curtailing the disease. The research extensively studies the SEQIAHR compartmental model of COVID-19 to provide insight into the dynamics of the disease by underlying tailored strategies designed to minimize the pandemic. We first studied the noncontrol model’s dynamic behaviour by calculating the reproduction number and examining the two nonnegative equilibria’ existence. The model utilizes the Castillo-Chavez method and Lyapunov function to investigate the global stability of the disease at the disease-free and endemic equilibrium. Sensitivity analysis was carried on to determine the impact of some parameters onR0. We further examined the COVID model to determine the type of bifurcation that it exhibits. To help contain the spread of the disease, we formulated a new SEQIAHR compartmental optimal control model with time-dependent controls: personal protection and vaccination of the susceptible individuals. We solved it by utilizing Pontryagin’s maximum principle after studying the dynamical behaviour of the noncontrol model. We solved the model numerically by considering different simulation controls’ pairing and examined their effectiveness.


---

## Body

## 1. Introduction

The unusual, life-threatening pneumonia condition affecting humanity remains the globe’s concern as governments struggle to defeat the pandemic. COVID-19, one of the most media campaigned viral diseases, emerged from Wuhan’s city in China in the last quarter of 2019. The epidemic has transcended the nations’ regions, and all sectors of the global economy have felt its adverse effects [1]. The epidemic has been the source of setbacks in businesses, disruption in academic calendars, and closure of production companies and public services [2]. The campaign against the disease has intensified due to the number of case counts and deaths recorded each day globally. The primary mode of transmissibility of the disease is airborne. Thus, individuals become infected by contact with the contaminated sneezes and droplets of the virus [3]. The documented number of COVID-19-related deaths and cases in the first quarter of the year 2020 caused governments to respond accordingly by promoting and observing self-protective protocols and restrictions directed by the World Health Organization. However, the protocols and restrictions were gradually eased when evidence of the disease waning was confirmed by governments [4]. Notably, the vast majority of the European countries are currently battling with the second wave of COVID-19 outbreaks after successfully minimizing COVID-19 disease in the early months of the year 2020 [5]. The second wave has caused the European countries’ governments to respond to the second spike by introducing new measures that may help manage the disease and prevent the further spread of the outbreak. These measures range from strict new restrictions, such as minilockdown and compulsory mask wearing, to relaxed regulations, allowing the opening of bars and pubs with a set time for closure. These recommended restrictions vary from country to country regarding the number of cases and infection rates [6].Since the detection of the virus in Wuhan, variants of COVID-19 mathematical models have been considered that have helped shape the pandemic, with inked preventive measures which could be adopted to flatten the curve, estimating the final epidemic size and prevent further infection (see [7–23]).Mathematical modelling of infectious diseases has been an experimenting station where the vast knowledge of epidemic outbreaks, transmission dynamics, intervention, and alternative courses of action for controlling diseases are discovered. Integrating mathematical computation into the modelling of infectious diseases has evolved in an unparalleled achievement regarding intervention strategies and has been the spine of all explorations on infectious diseases [24].In their paper, Hellewell et al. [25] assessed the effectiveness of a stochastic transmission model to control the new SARS-CoV-2 disease by utilizing the preventive measures of isolation and contact tracing. Qianying et al. [26], in their work, constructed SEIRNDC compartmental model of coronavirus-2 in Wuhan that examines the government’s responses towards the disease and the reaction from the people. The authors in [27] considered an eight compartmental nonlinear differential equation model of COVID-19 that uses time-dependent diagnosis and contact rates to refit their existing SEIR compartmental model [28] to a newly available data for better estimation of the basic reproduction number. In a related paper, Kang et al. [29] examined the COVID-19 epidemic dynamics in China. The authors applied Moran’s I spatial statistics in a test to ascertain the possibility of the disease’s spatial association. In [30], in their paper, the author adapted the growth model to study the spread pattern of MERS, COVID-19, and SARS by using the inhibition and growth rates to establish the model’s propagation. Benvenuto et al. [31] predicted the spread and trend of COVID-19 disease by using the ARIMA model. In a paper related to the current studies, Li et al. [32] investigated the COVID-19 transmission process using the official modelling data. The study examines the error between the model and the official data. In [7], in their work, the authors calibrated a COVID-19 epidemic model to deduce various characteristics such as age-dependent attack rates, length of incubation, generation periods, reproductive number, and growth rate of cumulative incidence. In [10], in their research, the authors analyzed the spread pattern of the COVID-19 outbreak by using a time series model. The authors in [33] considered variants of the SIR model with a parameter that factors into account the effects of social distancing. Fu et al. [34] applied Boltzmann-function-based regression analyses to estimate the number of SARS-CoV-2 confirmed cases in China. Shi et al. [35] considered a SEIR COVID-19 propagation model that assessed the variation in the length of the incubation period, weakness of the transmission ability of the incubation period, and the intervention of the government to detect and isolate the infected. In [36], the authors studied and predicted the pattern of the COVID-19 epidemic and estimated the various parameters involved in the model’s analyses. Liu et al. [37] considered a new QSEIR COVID-19 epidemic model that studied the effect of quarantine measures imposed in Mainland China during the peak of the outbreak. In [38], in their paper, a stochastic COVID-19 transmission model is examined that analyzed the intervention measures employed in Mainland China. Kissler et al. [39] applied a mathematical model to prove that the United States’ critical care capacity will not be sufficient to be maintained with a one-time intervention strategy. In [40], the authors fitted a coronavirus-2 age-structured model to data from six countries affected by the pandemic to determine the age gradient in observed cases. In [41], the authors studied a SEIR compartmental model of COVID-19, with modifications to account for the spreading of the epidemic in the latent stage and the effects of varying containment’s proportions. Acuna et al. [42] considered mathematical models which examine the issues related to the spread of SARS-CoV-2 in Mexico and explore the effect of behaviour changes needed to wane the infection. Rong et al. [43] studied a new dynamical model of SSqEI1I2HRV for SARS-CoV-2 to assess the effect of delay in detecting an infected person. In a related article, Giordano et al. [44], considered a nonlinear ordinary differential equation model of SIDARTHE compartments that predicts the propagation of COV-19 epidemic in Italy and helps to diffuse the misperceptions of epidemic spread and case fatality rate. In [45], the authors in their work explored the impact of community mask wearing on COVID-19 transmission dynamics in the United States by formulating a new mathematical model for the assessment.Optimal control models have engineered tailored strategies that have been paramount for minimizing and prevention of infections [46–56].Asamoah et al. [57] applied an optimal control theory to nonlinear ordinary differential equations of SEAIRV compartmental model of coronavirus transmission that analyzed the cost-effective strategy of all the proposed methods. In the paper by [58], the authors constructed and researched a compartmental SIRU model that brings further insight into the propagation of the COVID-19 disease. The SIRU compartmental model is then converted to an optimal control problem, and the numerical solutions are presented. The authors in [59] formulated a mathematical model for coronavirus-2 disease to estimate the model parameters by fitting it to actual data. The authors further performed optimal control analysis on the modified model. In [30], the author constructed a mathematical model that is segregated into compartments of SEIRV for the coronavirus-2 epidemic and converted it to an optimal control problem by characterizing several control strategies by applying the maximum principle. Sasmita et al. [60] considered SEI1I2RS nonlinear compartmental model of coronavirus-2 infection to predict the disease’s peak in Indonesia. The authors considered five time-dependent controls in constructing the optimal control model to deduce strategies critical for policymakers to curb the pandemic. The authors in [61] applied an optimal control analysis to a mathematical model of SARS-CoV-19 to help deduce many possible strategies for the control of the disease. In [62], the authors developed an ordinary differential equation model of SEIRW compartments that describes the COVID-19 disease’s dynamics. The authors made a SEIRW compartmental model modification to convert the dynamical model into a new SEIRW control model to assess the chosen strategies thoroughly.This research article presents a SEQIAHR compartmental model of COVID-19 to provide insight into the disease’s dynamics by utilizing tailored strategies to minimize the pandemic. The study is motivated by the available COVID-19 works and formulating a new SEQIAHR compartmental optimal control model that would add to the existing knowledge and help improve public health decision-making by providing scientific strategies to prevent the disease.The remaining work is organized as follows: Section2 formulates a deterministic SEQIAHR compartmental COVID-19 model. Section 3 studies the global stability of the model, sensitivity and bifurcation analysis. A new SEQIAHR compartmental optimal control model is constructed and analyzed in Section 4 with time-dependent control functions of personal protection and treatment of asymptomatic infected individuals. Finally, Section 5 discusses the simulated results of the models by using an iterative method of Runge-Kutta 4th-order method and Matlab.


## 2. The Model Formulation

This section formulates a compartmental SEQIAHR transmission model for COVID-19 disease to understand the dynamical behaviour of the disease and the strategy needed in curtailing it. Here, we modify the baseline model of [45] which is segregated into compartments of: susceptible St, exposed Et, symptomatic infectious It, asymptomatic infectious At, hospitalized Ht, recovered Rt, and cumulative deaths Dt, by considering recruitment into the population, Λ and natural death rate μ, and ignoring the cumulative death compartment. Additionally, we assume that COVID-19 has a high level of transmission, and therefore, the main way of preventing the spread of the disease is to adopt a quarantine approach of the exposed individuals. This is considered in the formulation of the model by first modifying the original model system (1) of [45], with the inclusion of quarantine compartment, by assuming that the exposed individuals are quarantined at rate τ. Further, the model assumes that the quarantine individuals are hospitalized at a rate of τ1. Also, we assume that the quarantined individuals may die naturally. We assumed a time-dependent population for the modified model. With these assumptions, the modified model equations are given as follows:
(1)ddtS=Λ−βSIN−βSηAN−μS,ddtE=βSIN+βSηAN−σE−τE−μE,ddtQ=τE−τ1+μQ,ddtI=ασE−ϕI−γII−μI,ddtA=1−ασE−γAA−μA,ddtH=ϕI+τ1Q−δH−γHH−μH,ddtR=γII+γAA+γHH−μR,withS≥0, E≥0, Q≥0, I≥0, A≥0, H≥0, and R≥0.


### 2.1. Analysis of Model: Positivity and Boundedness

Theorem 1.
The setSt,Et,Qt,It,At,Ht,Rt being the solution of the state Equation (1) with nonnegative parameters is positive with initial conditions given by the following:
(2)S0≥0,E0≥0,Q0≥0,I0≥0,A0≥0,H0≥0,R0≥0.Proof.
With the method illustrated by [57], theorem one can be proofed by adopting the same approach as underlaid below.
We redefineW=S,E,Q,I,A,H,RT, U0=β/NI, and U1=βη/NA, with T denoting the transposition. Hence, COVID-19 model Equation (1) in matrix form is represented as follows:
(3)dWdt=PW+B,where
(4)P=−U0+U1+μ000000U0+U1−σ+τ+μ000000τ−τ1+μ00000ασ0−ϕ+γI+μ00001−ασ00−γA+μ0000τ1ϕ0−δ+γH+μ0000γIγAγH−μ,B=Λ000000.In model Equation (1), rewriting the third equation into a first-order nonhomogenous differential equation gives
(5)dQdt+τ1+μQ=τE.Now, from (5), adopting the method of integrating factor and applying it provides
(6)Qt=e−τ1+μtQ0+τ∫0tEse−τ1+μs.ds.Similarly, mimicking the same approach, the fourth equation of model (1) gives
(7)It=e−ϕ+γI+μtI0+ασ∫0tEse−ϕ+γI+μs.ds.Now, as can be determined in the work of [57], it follows that by intuition, d/dtQ≥0, at t=0, and d/dtI≥0, at t=0, for E0=0. Hence, following the same steps, the same can be generalized for St, Et, At, Ht, and Rt, which ensures that the state variables stays positive in the entire time of the study. In addition, by inspection, it could be verified that the off-diagoanl entries of P are nonnegative and B≥0, which confirms the property Metzler matrix [57]. Hence, the state model Equation (1) is positively invariant in R+7.Theorem 2.
The nonlinear equation model (1) has solutions bounded within the invariant region, φ∈R7 given as
(8)φ=S,E,Q,I,A,H,R∈R+7,S+E+Q+I+A+H+R≤Λ−δH−μN.Proof.
Here, we denoteNt=S+E+Q+I+A+H+R. Then, the simplified nonlinear differential equation of (1) is given as
(9)ddtNt=Λ−δH−μN.
Then, from Equation (9), setting Z to be a solution of (9) gives a unique initial value problem
(10)ddtZ1t=Λ−μZ1tfort≥0,Z10=N0.
Solving Equation (54) produces
(11)Z1t=N0e−μt+Λμ1−e−μt.
Hence, from the comparison theorem [63], it can be seen that
(12)Nt=N0e−μt+Λμ1−e−μt.
Therefore, from (12), it follows that the possible solution set of the state variables S,E,Q,I,A,H,R is bounded and the state model Equation (1) is positively in variant in R+7, implying that any trajectory with initial condition in R+7 will remain in φ. This guarantees that model (1) is mathematically and epidemiologically feasible and well posed.


### 2.2. Disease-Free Equilibrium and Reproduction Ratio

The basic reproduction number threshold parameter is key in determining whether an invading infection from an infected individual will have the potency to be endemic or die out in a naive population. It is defined as a new infection from an infected individual who enters into a naive susceptible population. WhenR=0, H=0, A=0, I=0, Q=0, and E=0, then the COVID-19 model (1) has a disease-free equilibrium (DFE), determined by equating the right-hand side of the equations in model (1) to zero and solve for the state variables. Hence,
(13)E0=S0,E0,Q0,I0,A0,H0,R0=Λμ,0,0,0,0,0,0.The basic reproduction number of the model system (1) at E0=S0,E0,Q0,I0,A0,H0,R0 is derived from the method studied in [64] by Diekmann et al. Based on [64], we derive matrices F and V as follows:
(14)FJ=00ββη000000000000,V=σ+τ+μ000−ττ1+μ00−ασ0ϕ+γI+μ0−1−ασ00γA+μ.Hence, the basic reproduction number for the COVID-19 model (1) is obtained by the spectral radius of ρFV−1 as follows:
(15)R0=ασβσ+τ+μϕ+γI+μ+1−αησβσ+τ+μγA+μ.


### 2.3. Existence of Endemic Equilibrium

WhenR≠0, H≠0, A≠0, I≠0, Q≠0, E≠0, and S≠0, then it follows that the endemic equilibrium exists. This unique equilibrium point for the COVID-19 model system (1) given by E1∗=S∗,E∗,Q∗,I∗,A∗,H∗,R∗ is determined as follows:
(16)S∗=1R0,E∗=μNR0−1ϕ+γI+μγA+μNβσαγA+μ+η1−αϕ+γI+μ,Q∗=μτ1NR0−1ϕ+γI+μγA+μNβσξ+μαγA+μ+η1−αϕ+γI+μ,I∗=ασμNR0−1γA+μNβσαγA+μ+η1−αϕ+γI+μ,A∗=1−ασμNR0−1ϕ+γI+μNβσαγA+μ+η1−αϕ+γI+μ,H∗=NR0−1γA+μμNϕασ+τ1ϕ+γI+μβσαγA+μ+η1−αϕ+γI+μ,R∗=σμNR0−1αγIγA+μ+γA1−αϕ+γI+μ+αϕηγHγA+μβσαγA+μ+η1−αϕ+γI+μ.


## 2.1. Analysis of Model: Positivity and Boundedness

Theorem 1.
The setSt,Et,Qt,It,At,Ht,Rt being the solution of the state Equation (1) with nonnegative parameters is positive with initial conditions given by the following:
(2)S0≥0,E0≥0,Q0≥0,I0≥0,A0≥0,H0≥0,R0≥0.Proof.
With the method illustrated by [57], theorem one can be proofed by adopting the same approach as underlaid below.
We redefineW=S,E,Q,I,A,H,RT, U0=β/NI, and U1=βη/NA, with T denoting the transposition. Hence, COVID-19 model Equation (1) in matrix form is represented as follows:
(3)dWdt=PW+B,where
(4)P=−U0+U1+μ000000U0+U1−σ+τ+μ000000τ−τ1+μ00000ασ0−ϕ+γI+μ00001−ασ00−γA+μ0000τ1ϕ0−δ+γH+μ0000γIγAγH−μ,B=Λ000000.In model Equation (1), rewriting the third equation into a first-order nonhomogenous differential equation gives
(5)dQdt+τ1+μQ=τE.Now, from (5), adopting the method of integrating factor and applying it provides
(6)Qt=e−τ1+μtQ0+τ∫0tEse−τ1+μs.ds.Similarly, mimicking the same approach, the fourth equation of model (1) gives
(7)It=e−ϕ+γI+μtI0+ασ∫0tEse−ϕ+γI+μs.ds.Now, as can be determined in the work of [57], it follows that by intuition, d/dtQ≥0, at t=0, and d/dtI≥0, at t=0, for E0=0. Hence, following the same steps, the same can be generalized for St, Et, At, Ht, and Rt, which ensures that the state variables stays positive in the entire time of the study. In addition, by inspection, it could be verified that the off-diagoanl entries of P are nonnegative and B≥0, which confirms the property Metzler matrix [57]. Hence, the state model Equation (1) is positively invariant in R+7.Theorem 2.
The nonlinear equation model (1) has solutions bounded within the invariant region, φ∈R7 given as
(8)φ=S,E,Q,I,A,H,R∈R+7,S+E+Q+I+A+H+R≤Λ−δH−μN.Proof.
Here, we denoteNt=S+E+Q+I+A+H+R. Then, the simplified nonlinear differential equation of (1) is given as
(9)ddtNt=Λ−δH−μN.
Then, from Equation (9), setting Z to be a solution of (9) gives a unique initial value problem
(10)ddtZ1t=Λ−μZ1tfort≥0,Z10=N0.
Solving Equation (54) produces
(11)Z1t=N0e−μt+Λμ1−e−μt.
Hence, from the comparison theorem [63], it can be seen that
(12)Nt=N0e−μt+Λμ1−e−μt.
Therefore, from (12), it follows that the possible solution set of the state variables S,E,Q,I,A,H,R is bounded and the state model Equation (1) is positively in variant in R+7, implying that any trajectory with initial condition in R+7 will remain in φ. This guarantees that model (1) is mathematically and epidemiologically feasible and well posed.


## 2.2. Disease-Free Equilibrium and Reproduction Ratio

The basic reproduction number threshold parameter is key in determining whether an invading infection from an infected individual will have the potency to be endemic or die out in a naive population. It is defined as a new infection from an infected individual who enters into a naive susceptible population. WhenR=0, H=0, A=0, I=0, Q=0, and E=0, then the COVID-19 model (1) has a disease-free equilibrium (DFE), determined by equating the right-hand side of the equations in model (1) to zero and solve for the state variables. Hence,
(13)E0=S0,E0,Q0,I0,A0,H0,R0=Λμ,0,0,0,0,0,0.The basic reproduction number of the model system (1) at E0=S0,E0,Q0,I0,A0,H0,R0 is derived from the method studied in [64] by Diekmann et al. Based on [64], we derive matrices F and V as follows:
(14)FJ=00ββη000000000000,V=σ+τ+μ000−ττ1+μ00−ασ0ϕ+γI+μ0−1−ασ00γA+μ.Hence, the basic reproduction number for the COVID-19 model (1) is obtained by the spectral radius of ρFV−1 as follows:
(15)R0=ασβσ+τ+μϕ+γI+μ+1−αησβσ+τ+μγA+μ.


## 2.3. Existence of Endemic Equilibrium

WhenR≠0, H≠0, A≠0, I≠0, Q≠0, E≠0, and S≠0, then it follows that the endemic equilibrium exists. This unique equilibrium point for the COVID-19 model system (1) given by E1∗=S∗,E∗,Q∗,I∗,A∗,H∗,R∗ is determined as follows:
(16)S∗=1R0,E∗=μNR0−1ϕ+γI+μγA+μNβσαγA+μ+η1−αϕ+γI+μ,Q∗=μτ1NR0−1ϕ+γI+μγA+μNβσξ+μαγA+μ+η1−αϕ+γI+μ,I∗=ασμNR0−1γA+μNβσαγA+μ+η1−αϕ+γI+μ,A∗=1−ασμNR0−1ϕ+γI+μNβσαγA+μ+η1−αϕ+γI+μ,H∗=NR0−1γA+μμNϕασ+τ1ϕ+γI+μβσαγA+μ+η1−αϕ+γI+μ,R∗=σμNR0−1αγIγA+μ+γA1−αϕ+γI+μ+αϕηγHγA+μβσαγA+μ+η1−αϕ+γI+μ.


## 3. Disease-Free Equilibrium (DFE) and Its Stability

In this section, the global stability analysis of the COVID-19 model (1) at the disease-free equilibrium is studied. The Castillo-Chavez method [65] would be used to prove that model (1) is globally asymptotically stable at the disease-free equilibrium. Thus, it follows that in considering the method of Castillo-Chavez, the COVID-19 model (1) can be transformed as follows:
(17)dp1dt=y1p1,p2,dp2dt=y2p1,p2,y2p1,0=0,wherep1 denotes the uninfected population; thus, p1=S,R, and p2 represents the infected, with, p2=E,Q,I,A,H. The disease-free equilibrium point of (1) is given by U=p10,0.The pointy10,0 is a globally stable asymptotically stable equilibrium for the model (1) provided R0<1, and the below criteria are satisfied.D1. Given dp1/dt=y1p1,0, p10 is globally asymptotically stable.D2. yp1,p2=Zp2−y^2p1,p2, where y^2p1,p2≥0 for p1,p2∈ζu.What happens next is if the model Equation (1) meets the above conditions; then, the following theorem holds.Theorem 3.
The pointU=p10,0 is globally asymptotically stable equilibrium given that R0<1 and the conditions D1 and D2 are satisfied.Proof.
Concerning the model Equation (1), we derive y1p1,p2 and y2p1,p2 as
(18)y1p1,p2=Λ−βSIN−βSηAN−μSγII+γAA+γHH−μR,y2p1,p2=βSIN+βSηAN−σE−τE−μEτE−τ1+μQασE−ϕI−γII−μI1−ασE−γAA−μAϕI+τ1Q−δH−γHH−μH.
It follows that withS=S0, I=I0, A=A0, H=H0, and R=R0, then y1p1,0 becomes
(19)y1p1,p2=Λ−βS0I0N−βS0ηA0N−μS0γII0+γAA0+γHH0−μR0.
What follows next is that from (19), we notice that as t⟶∞, p1=p10. Hence, p1=p10 is globally asymptotically stable, which verifies the first condition.
Now, in determining whether condition two would be satisfied, we utilizeYp1,p2=Zp2−y^2p1,p2. And we get
(20)A=−a110βNS0βηNS00τ−a22000ασ0−a33001−ασ00−a4400τ1ϕ0−a55EQIAH−y^2p1,p20000,
whereA=Dy2−W^2y1,y2, a11=σ+τ+μ, a22=τ1+μ, a33=ϕ+γI+μ, a44=γA+μ, a55=δ+γH+μ, and matrix D given by
(21)Z=−a110βNS0βηNS00τ−a22000ασ0−a33001−ασ00−a4400τ1ϕ0−a55,
with(22)J=βS0−SINβηS0−SIN00000000,
whereJ=y^2p1,p2. It can be ascertained from model Equation (1) that the total population is bounded by S0. Therefore, it follows that βI/NS≤βI/NS0, and ηβI/NS≤ηβI/NS0 which implies y^2p1,p2 is positive definite. Further, matrix Z is evidently an M-matrix, with the off-diagonal entries positive. Hence, condition two is satisfied which proves the global asymptotic stability of U.


### 3.1. Endemic Equilibrium (EE) and Its Stability

This subsection presents the global stability analysis of the COVID-19 model (1) at the endemic equilibrium by applying a Lyapunov function theory [66] for the global stability analysis. The results are given as follows:Theorem 4.
The unique endemic equilibriumE1∗ for the COVID-19 model (1) is globally asymptotically stable in R+7 whenever R0>1.Proof.
The Lyapunov functionL=m1S+m2E+m3Q+m4I+m5A+m6H+m7R, where mi for i=1,2,3⋯.7 are constants to be chosen in the course of the proof are defined.
The derivative ofL along the solution of (1) is given by
(23)dLdt=m1dSdt+m2dEdt+m3dqdt+m4dIdt+m5dAdt+m6dHdt+m7dRdt,dLdt=m1Λ−βI+ηASN−μS+m2βI+ηASN−σE−τE−μE+m3τE−τ1+μQ+m4ασE−ϕI−γII−μI+m51−ασE−γAA−μA+m6ϕI+τ1Q−δH−γHH−μH+m7γII+γAA+γHH−μR,dLdt=m1Λ−μS+m2−m1βI+ηASN+m5−m2σE+m3−m2τE−m2μE+m6−m3τ1Q−m3μQ+m4−m5ασE+m6−m4ϕI+m7−m4γII−m4μI+m7−m5γAA−m5μA−m6δH−m6μH+m7−m6γHH−m7μR.
Choosingm1, m2, m3, m4, m5, m6, and m7 such that m1=m2=m3=m4=m5=m6=m7 and Λ−μS=0 gives
(24)−m2μE−m3μQ−m4μI−m5μA−m6δH−m6μH−m7μR.
It follows thatL is positive definite, and dL/dt is negative definite. Therefore, the function L is a Lyapunov function for model system (1), and by Lyapunov asymptotic stability theorem [67], the endemic equilibrium E1∗ is globally asymptotically stable.


### 3.2. Sensitivity Analysis

Sensitivity analysis is one of the essential subjects that has been explored by many researchers and is of great importance to epidemiological modelling. Sensitivity analysis study assists us in ascertaining parameters that impact theR0 and allows epidemiologists to improve the design of the control strategies. The results from sensitivity index computation indicate the effect of involving parameters that contribute to the spreading of the epidemic and inform us of the relative change of Ro and other parameters.Definition 5.
For a given parameterα, the normalized forward sensitivity index of R0 is computed using the formula discussed in [18, 68], as
(25)ςαR0=∂R0∂ααR0.
Applying this formula for the parametersβ and α gives
(26)∂R0∂ββR0=ασγA+μ+1−αησϕ+γI+μσ+τ+μϕ+γI+μγA+μ.βσ+τ+μϕ+γI+μγA+μασβγA+μ+1−αησβϕ+γI+μ=1,∂R0∂ααR0=1−ασβσ+τ+μγA+μ.ησ+τ+μϕ+γI+μγA+μασβγA+μ+1−αησβϕ+γI+μ=1−ασηασβγA+μ+1−αησβϕ+γI+μ.Mimicking the above method for the remaining parameters and evaluating the results with the parameter values of Table1 provide the sensitivity indices of R0 parameters presented in Table 2.Table 1
COVID-19 model 1 equation parameters.
ParameterDescriptionValueReferenceβInfection contact rate(1.5)/day[26]σTransition from exposed to infectious(1/14)/day[26]ηInfectiousness factor for asymptomatic carriers(0.6)[26]αFraction of infections that become symptomatic(0.15)[26]ϕHospitalization rate(0.02)[26]γAAsymptomatic (recovery rate)(1/14)/day[26]γISymptomatic (recovery rate)(1/30)/day[26]γHHospitalized (recovery rate)(1/14)/day[66]δDeath rate (hospitalized)0.01[45]ΛRecruitment rate50AssumedτQuarantine rate0.012[60]τ1Hospitalized rate0.06[26]μNatural death rate0.000042578[4]Table 2
Parameters forR0 and their sensitivity index for model (1).
ParameterSensitivity indexβ+1σ-1.3094η+0.7174α+0.1560ϕ-0.1066γA-0.0510γI-0.1756τ-0.1445μ-0.00010As noted from Table2, β, η, and α are the parameters with positive indices contributing to the spreading of the epidemic. The positive parameters contribute to the spreading of the outbreak since they increase the R0. On the other hand, the parameters with a negative index contribute to controlling the disease since they have reduced the R0. Further, as noted, the parameter β has a sensitivity index of +1, which implies increasing or decreasing β by a specific percentage increases or reduces R0 by the same percentage.


### 3.3. Bifurcation Analysis

In an attempt to ascertain whether model system (1) exhibits backward bifurcation or not, we analyze model system (1) with the theory of centre manifold as credited to Castillo-Chavez and Song in ([69] see Theorem 6).Hence, we use the approach of ([69, 70]) to determine the criteria on which the parameter values of model system (1) cause a backward or forward bifurcation to occur. We consider the system below:
(27)ddtx=fx,∈1,wheref is continously differentiable at least twice in x and ∈1 is the bifurcation parameters.Equationsa and b are denoted by as follows:
(28)a=∑k,i,j=1nvkwiwj∂2fk0,0∂xi∂xjand
(29)b=∑k,i,j=1nvkwi∂2fk0,0∂xi∂∈1are the determinants of the existence of bifurcation in model system (1). Thus, when a>0 and b>0, then backward bifurcation exists and occurs at ∈1=0. When a<0 and b>0, forward bifurcation exists and occurs at ∈1=0. Now, considering β as a bifurcation parameter, then R0 be equivalent to
(30)β=β∗=σ+τ+μϕ+γI+μγA+μασγA+μ+1−αησσ+τ+μ.We alter the initial design of the state model (1) for easy computations as follows; S=x1, E=x2, Q=x3, I=x4, A=x5, H=x6, and R=x7, such that N=x1+x2+x4+x5+x7. Also, X=x1,x2,⋯,x7T and f=f1,f2,⋯,f7T are vector notations of the model system (1) and can be rewritten in the form
(31)ddtx=fx,β∗,where
(32)∑i=17fi=ddtxi=ddtx1=Λ−βx1x4∑i=1nxi−βx1ηx5∑i=1nxi−μx1,ddtx2=βx1x4∑i=1nxi+βx1ηx5∑i=1nxi−σE−τx2−μx2,ddtx3=τx2−τ1+μx3,ddtx4=ασx2−ϕx4−γIx4−μx4,ddtx5=1−ασx2−γAx5−μx5,ddtx6=τ1x3+ϕx4−δx6−γHx6−μx6,ddtx7=γIx4+γAx5+γHx6−μx7.Now, the Jacobian matrix of model system (1) was evaluated at E0, when β=β∗ is provided as
(33)JE0β=β∗=−μ00−β∗−β∗η000−σ+μ0−β∗−β∗η000−τ−τ1+μ00000ασ0−ϕ+γI+μ00001−ασ00−γA+μ0000τ1ϕ0−δ+γH+μ0000γIγAγH−μ.The Jacobian matrixJE0β=β∗ has a right eigenvector corressponding to a simple zero eigenvalue given by W=w1,w2,⋯,w7T and a left eigenvector corressponding to a simple eigenvalue given by V=v1,v2,⋯,v7T. Deducing the eigenvectors of W and V gives
(34)w1=−ασβ∗w2μϕ+γI+μ−1−ασβ∗w2γA+μ,w3=−τw2τ1+μ,w2>0,w4=ασw2ϕ+γI+μ,w5=1−ασw2γA+μ,w6=−ττ1w2δ+γH+μτ1+μ+ϕασw2ϕ+γI+μ,w7=ασγIw2μϕ+γI+μ+1−ασγAw2μγA+μ+ϕασγHw2μϕ+γI+μ−ττ1γHw2μδ+γH+μτ1+μ,(35)v1=v7=0,v2>0,v3=0,v6=0,v4=β∗v2ϕ+γI+μ,v5=β∗v2γA+μ.The derivation of the nonzero partial derivatives is given by(36)∂2f2∂x1∂x4=−αv2β∗2σ2w22αγA+μ+1−αϕ+γI+μμμϕ+γI+μ2γA+μ,∂2f2∂x1∂x5=−v2β∗2ησ2w221−ααγA+μ+1−αϕ+γI+μμμϕ+γI+μ2γA+μ,∂2f2∂x4∂β∗=v2ασΛw2μϕ+γI+μ,∂2f2∂x5∂β∗=Λv21−ασw2μγA+μ.Hence, we obtain(37)a=−v2β∗2ησ2w221−ααγA+μ+1−αϕ+γI+μμαγA+μ+1−αϕ+γI+μμμϕ+γI+μ2γA+μ,b=v2σΛw2αγA+μ+1−αϕ+γI+μμμϕ+γI+μγA+μ.The coefficientb is positive as always. According to Theorem 6 of Castillo-Chavez and Song [69], the sign of a determines the local dynamics around the disease-free equilibrium for β=β∗.


## 3.1. Endemic Equilibrium (EE) and Its Stability

This subsection presents the global stability analysis of the COVID-19 model (1) at the endemic equilibrium by applying a Lyapunov function theory [66] for the global stability analysis. The results are given as follows:Theorem 4.
The unique endemic equilibriumE1∗ for the COVID-19 model (1) is globally asymptotically stable in R+7 whenever R0>1.Proof.
The Lyapunov functionL=m1S+m2E+m3Q+m4I+m5A+m6H+m7R, where mi for i=1,2,3⋯.7 are constants to be chosen in the course of the proof are defined.
The derivative ofL along the solution of (1) is given by
(23)dLdt=m1dSdt+m2dEdt+m3dqdt+m4dIdt+m5dAdt+m6dHdt+m7dRdt,dLdt=m1Λ−βI+ηASN−μS+m2βI+ηASN−σE−τE−μE+m3τE−τ1+μQ+m4ασE−ϕI−γII−μI+m51−ασE−γAA−μA+m6ϕI+τ1Q−δH−γHH−μH+m7γII+γAA+γHH−μR,dLdt=m1Λ−μS+m2−m1βI+ηASN+m5−m2σE+m3−m2τE−m2μE+m6−m3τ1Q−m3μQ+m4−m5ασE+m6−m4ϕI+m7−m4γII−m4μI+m7−m5γAA−m5μA−m6δH−m6μH+m7−m6γHH−m7μR.
Choosingm1, m2, m3, m4, m5, m6, and m7 such that m1=m2=m3=m4=m5=m6=m7 and Λ−μS=0 gives
(24)−m2μE−m3μQ−m4μI−m5μA−m6δH−m6μH−m7μR.
It follows thatL is positive definite, and dL/dt is negative definite. Therefore, the function L is a Lyapunov function for model system (1), and by Lyapunov asymptotic stability theorem [67], the endemic equilibrium E1∗ is globally asymptotically stable.


## 3.2. Sensitivity Analysis

Sensitivity analysis is one of the essential subjects that has been explored by many researchers and is of great importance to epidemiological modelling. Sensitivity analysis study assists us in ascertaining parameters that impact theR0 and allows epidemiologists to improve the design of the control strategies. The results from sensitivity index computation indicate the effect of involving parameters that contribute to the spreading of the epidemic and inform us of the relative change of Ro and other parameters.Definition 5.
For a given parameterα, the normalized forward sensitivity index of R0 is computed using the formula discussed in [18, 68], as
(25)ςαR0=∂R0∂ααR0.
Applying this formula for the parametersβ and α gives
(26)∂R0∂ββR0=ασγA+μ+1−αησϕ+γI+μσ+τ+μϕ+γI+μγA+μ.βσ+τ+μϕ+γI+μγA+μασβγA+μ+1−αησβϕ+γI+μ=1,∂R0∂ααR0=1−ασβσ+τ+μγA+μ.ησ+τ+μϕ+γI+μγA+μασβγA+μ+1−αησβϕ+γI+μ=1−ασηασβγA+μ+1−αησβϕ+γI+μ.Mimicking the above method for the remaining parameters and evaluating the results with the parameter values of Table1 provide the sensitivity indices of R0 parameters presented in Table 2.Table 1
COVID-19 model 1 equation parameters.
ParameterDescriptionValueReferenceβInfection contact rate(1.5)/day[26]σTransition from exposed to infectious(1/14)/day[26]ηInfectiousness factor for asymptomatic carriers(0.6)[26]αFraction of infections that become symptomatic(0.15)[26]ϕHospitalization rate(0.02)[26]γAAsymptomatic (recovery rate)(1/14)/day[26]γISymptomatic (recovery rate)(1/30)/day[26]γHHospitalized (recovery rate)(1/14)/day[66]δDeath rate (hospitalized)0.01[45]ΛRecruitment rate50AssumedτQuarantine rate0.012[60]τ1Hospitalized rate0.06[26]μNatural death rate0.000042578[4]Table 2
Parameters forR0 and their sensitivity index for model (1).
ParameterSensitivity indexβ+1σ-1.3094η+0.7174α+0.1560ϕ-0.1066γA-0.0510γI-0.1756τ-0.1445μ-0.00010As noted from Table2, β, η, and α are the parameters with positive indices contributing to the spreading of the epidemic. The positive parameters contribute to the spreading of the outbreak since they increase the R0. On the other hand, the parameters with a negative index contribute to controlling the disease since they have reduced the R0. Further, as noted, the parameter β has a sensitivity index of +1, which implies increasing or decreasing β by a specific percentage increases or reduces R0 by the same percentage.


## 3.3. Bifurcation Analysis

In an attempt to ascertain whether model system (1) exhibits backward bifurcation or not, we analyze model system (1) with the theory of centre manifold as credited to Castillo-Chavez and Song in ([69] see Theorem 6).Hence, we use the approach of ([69, 70]) to determine the criteria on which the parameter values of model system (1) cause a backward or forward bifurcation to occur. We consider the system below:
(27)ddtx=fx,∈1,wheref is continously differentiable at least twice in x and ∈1 is the bifurcation parameters.Equationsa and b are denoted by as follows:
(28)a=∑k,i,j=1nvkwiwj∂2fk0,0∂xi∂xjand
(29)b=∑k,i,j=1nvkwi∂2fk0,0∂xi∂∈1are the determinants of the existence of bifurcation in model system (1). Thus, when a>0 and b>0, then backward bifurcation exists and occurs at ∈1=0. When a<0 and b>0, forward bifurcation exists and occurs at ∈1=0. Now, considering β as a bifurcation parameter, then R0 be equivalent to
(30)β=β∗=σ+τ+μϕ+γI+μγA+μασγA+μ+1−αησσ+τ+μ.We alter the initial design of the state model (1) for easy computations as follows; S=x1, E=x2, Q=x3, I=x4, A=x5, H=x6, and R=x7, such that N=x1+x2+x4+x5+x7. Also, X=x1,x2,⋯,x7T and f=f1,f2,⋯,f7T are vector notations of the model system (1) and can be rewritten in the form
(31)ddtx=fx,β∗,where
(32)∑i=17fi=ddtxi=ddtx1=Λ−βx1x4∑i=1nxi−βx1ηx5∑i=1nxi−μx1,ddtx2=βx1x4∑i=1nxi+βx1ηx5∑i=1nxi−σE−τx2−μx2,ddtx3=τx2−τ1+μx3,ddtx4=ασx2−ϕx4−γIx4−μx4,ddtx5=1−ασx2−γAx5−μx5,ddtx6=τ1x3+ϕx4−δx6−γHx6−μx6,ddtx7=γIx4+γAx5+γHx6−μx7.Now, the Jacobian matrix of model system (1) was evaluated at E0, when β=β∗ is provided as
(33)JE0β=β∗=−μ00−β∗−β∗η000−σ+μ0−β∗−β∗η000−τ−τ1+μ00000ασ0−ϕ+γI+μ00001−ασ00−γA+μ0000τ1ϕ0−δ+γH+μ0000γIγAγH−μ.The Jacobian matrixJE0β=β∗ has a right eigenvector corressponding to a simple zero eigenvalue given by W=w1,w2,⋯,w7T and a left eigenvector corressponding to a simple eigenvalue given by V=v1,v2,⋯,v7T. Deducing the eigenvectors of W and V gives
(34)w1=−ασβ∗w2μϕ+γI+μ−1−ασβ∗w2γA+μ,w3=−τw2τ1+μ,w2>0,w4=ασw2ϕ+γI+μ,w5=1−ασw2γA+μ,w6=−ττ1w2δ+γH+μτ1+μ+ϕασw2ϕ+γI+μ,w7=ασγIw2μϕ+γI+μ+1−ασγAw2μγA+μ+ϕασγHw2μϕ+γI+μ−ττ1γHw2μδ+γH+μτ1+μ,(35)v1=v7=0,v2>0,v3=0,v6=0,v4=β∗v2ϕ+γI+μ,v5=β∗v2γA+μ.The derivation of the nonzero partial derivatives is given by(36)∂2f2∂x1∂x4=−αv2β∗2σ2w22αγA+μ+1−αϕ+γI+μμμϕ+γI+μ2γA+μ,∂2f2∂x1∂x5=−v2β∗2ησ2w221−ααγA+μ+1−αϕ+γI+μμμϕ+γI+μ2γA+μ,∂2f2∂x4∂β∗=v2ασΛw2μϕ+γI+μ,∂2f2∂x5∂β∗=Λv21−ασw2μγA+μ.Hence, we obtain(37)a=−v2β∗2ησ2w221−ααγA+μ+1−αϕ+γI+μμαγA+μ+1−αϕ+γI+μμμϕ+γI+μ2γA+μ,b=v2σΛw2αγA+μ+1−αϕ+γI+μμμϕ+γI+μγA+μ.The coefficientb is positive as always. According to Theorem 6 of Castillo-Chavez and Song [69], the sign of a determines the local dynamics around the disease-free equilibrium for β=β∗.


## 4. COVID-19 Optimal Control Model

In formulating the optimal control model, we restructure the compartmental model (1) into an optimal control model with admissible controls that are considered to be continuous in time. The controls that are identified for the new structured control model are defined as follows:F1. The personal protection control rate varies with time and is given by n1.F2. The vaccination control rate varies with time and is denoted as n2.F3. The considered time is given by t∈0,T, where T is the final time and relatively short.Hence, the new restructured system for Equation (1) is provided below.
(38)ddtS=Λ−1−n1βSIN−1−n1βSηAN−μS−n2S,ddtE=1−n1βSIN+1−n1βSηAN−σE−τE−μE,ddtQ=τE−τ1+μQ,ddtI=ασE−ϕI−γII−μI,ddtA=1−ασE−γAA−μA,ddtH=ϕI+τ1Q−δH−γHH−μH,ddtR=γII+γAA+γHH+n2S−μR,withS≥0, E≥0, Q≥0, I≥0, A≥0, H≥0, and R≥0.We usher in a measurable control set:(39)C≔n=n1,n2njtisLebesguemeasurable,0≤nit≤1,t∈0,tfforj=1,2..The target of the considered control strategy is to(a)
Lower the COVID exposed, asymptomatic, and symptomatic infectious persons(b)
Make intervention cost small as possibleIn achieving the intended goals, we design an objective functional below as in ([57, 71]).
(40)J=∫0tfν1E+ν2I+ν3A+12h1n12+12h2n22.dt.The constantsν1, ν2, and ν3 are weight related to exposed, symptomatic infectious, and asymptomatic infectious individuals, respectively. Additionally, the weights h1 and h2 are positive and in association with time-dependent control functions n1,n2, respectively.The main objective of the control mode is to identify an optimal control pairn∗=n1∗,n2∗ that makes
(41)Jn∗=minNJn1,n2.The fundamental concept of the optimal control problem requires that we verify the existence and uniqueness of the optimal controls to characterize them.


### 4.1. Existence of Optimal Controls

As noted in [72], the existence result of Fleming and Richel is considered to show the existence of optimal control duple that minimizes (40) subject to the system (38).Theorem 6.
An optimal control duplen∗ exists that minimizes the objective functional (40) subject to the system (38), given that the below properties are met.
(a)
The set of control is convex and closed(b)
The system (38) is bounded by a linear function in both the state and control variable(c)
The objective functional (40) integrand is convex with respect to the control(d)
There exist constantsb1, b2≥0, and b3≥1 that make the objective functional (40) integrand bounded by the below quantity(42)b1∑i=12ni2b3/2−b2.Proof.
(a)
It is sufficient to writen=n1×n2 by definition of (39). Hence, n=n1×n2 is bounded and convex ∀t∈0,T. In addition, we choose k,d∈N, so that k=k1,k2 and d=d1,d2. Then, ∀v∈0,1, and we have vki+1−vdi∈N, satisfying the convexity property of the control set(b)
We denote the right side of system (38) and the associated solution by f and θ; then,(43)ψ=−1−n1βS1I1N−S2I2N−1−n1βηS1A1N−S2A2N−μS1−S2−n2S1−S21−n1βS1I1N−S2I2N+1−n1βηS1A1N−S2A2N−σE1−E2−τE1−E2−μE1−E2τE1−E2−τ1+μQ1−Q2ασE1−E2−ϕI1−I2−γII1−I2−μI1−I21−ασE1−E2−γAA1−A2−μA1−A2ϕI1−I2+τ1Q1−Q2−δH1−H2−γHH1−H2−μH1−H2γII1−I2+γAA1−A2+n2S1−S2−μR1−R2,where ψ=ft,θ1,n−ft,θ2,n(44)≤2β1−n1S1I1N−S2I2N+2βη1−n1S1A1N−S2A2N+2τE1−E2+2σE1−E2+2ασE1−E2+2τ1Q1−Q2+2γAA1−A2+2γHH1−H2+2n2S1−S2+2ϕI1−I2+μS2−S1+μE2−E1+μQ2−Q1+μI2−I1+μA2−A1+μH2−H1+μR2−R1,≤2β1−n1I1S1−S2+S2I1−I2+2βη1−n1A1S1−S2+S2A1−A2+2τ+2σ+2α+2ασ+μE1−E2+2τ1+μQ1−Q2+2γA+μA1−A2+2γH+μH1−H2+2n2+μS1−S2+2ϕ+μI1−I2+μR2−R1,≤S1−S22β1−n1I1+2βη1−n1A1+2β1−n1S2I1−I2+2βη1−n1S2A1−A2+2τ+2σ+2α+2ασ+μE1−E2+2τ1+μQ1−Q2+2γA+μA1−A2+2γH+μH1−H2+2n2+μS1−S2+2ϕ+μI1−I2+μR2−R1,≤2n2+μS1−S2+2τ+2σ+2α+2ασ+μE1−E2+2τ1+μQ1−Q2+2βΛμ1−n1+2ϕ+μI1−I2+2βηΛμ1−n1+2γA+μA1−A2+2γH+μH1−H2+μR2−R1,≤G1S1−S2+G2E1−E2+G3Q1−Q2+G4I1−I2+G5A1−A2+G6H1−H2+G7R1−R2,
whereG1=2n2+μ, G2=2τ+2σ+2ασ+μ, G3=2τ1+μ, G4=2βΛ/μ1−n1+2ϕ+μ, G5=2βηΛ/μ1−n1+2γA+μ, G6=2γA+μ, G7=μ, and G=maxG1,G2,G3,G4,G5,G6,G7.
Hence,f is uniformly Lipschitz continuous.
(c)
The Langrangian defined asLt,z,n is the integrand of the objective functional (40). Thus, we rewrite Lt,z,n in the form(45)Lt,z,n=r1t,z+r2t,z,with
r1t,z=ν1E+ν2I+ν3A and r2t,z=1/2∑j=12hjnj. The convexity of r2t,z, which is a linear combination of the control function 1/2∑j=12hjnj, needs to be proved. We prove the convexity by letting c:0,12⟶R be hn=1/2n2. Then, ∀u1,u2∈0,12 and ρ∈0,1. Hence, it follows that the below inequality holds.
(46)ρcu1+1−ρcu2≥cρu1+1−ρu2.
This confirms the convexity of the Langrangian with respect to the control.(d)
When observed from (45), we easily see that Lt,z,n≥r1t,z. Hence, we conclude that(47)Lt,z,n≥12h1n12+12h2n22≥b1∑i=12ni2b3/2−b2,with b1=1/2minh1,h2, b2>0, and b3=2. This completes the proof.


### 4.2. Characterization of Optimal Controls

Pontryagin’s maximum principle has been the wheel on which the necessary condition for the COVID (38) duple control needs to meet. With the principle, we convert the COVID (38) and the objective functional (40) into a problem of minimizing the Hamiltonian Hf with respect to the controls njt,j=1,2. Hence, the Hamiltonian Hf is given by
(48)Hf=ν1E+ν2I+ν3A+12h1n12+12h2n22+ζ1Λ−1−n1βSIN−1−n1βSηAN−μS−n2S+ζ21−n1βSIN+1−n1βSηAN−σE−τE−μE+ζ3τE−τ1+μQ+ζ4ασE−ϕI−γII−μI+ζ51−ασE−γAA−μA+ζ6ϕI+τ1Q−δH−γHH−μH+ζ7γII+γAA+γHH+n2S−μR.Theorem 7.
With the optimal control duplen1∗,n2∗ satisfying the condition (41), there exist adjoint variables ζi satisfying the adjoint system below.
(49)dζ1dt=ζ1−ζ21−n1βIS+E+I+A+R+SN2+ζ1−ζ21−n1βηAS+E+I+A+R+SN2+ζ1−ζ7n2+μζ1,dζ2dt=−ν1+ζ2−ζ1βSIN2+ζ2−ζ11−n1βηSAN2+ζ2−ζ5σ+ζ2−ζ3τ+ζ5−ζ4ασ+μζ2,dζ3dt=ζ3−ζ6τ1+μζ3,dζ4dt=−ν2+ζ1−ζ21−n1βSS+E+I+A+R+IN2+ζ2−ζ11−n1βηSAN2+ζ4−ζ6ϕ+ζ4−ζ7γI+μζ4,dζ5dt=−ν3+ζ1−ζ21−n1βηAS+E+I+A+R+AN2+ζ2−ζ1βSIN2+ζ5−ζ7γA+μζ5,dζ6dt=ζ6−ζ7γH+δζ6+μζ6,dζ7dt=ζ2−ζ1βSIN2+ζ2−ζ1βηSAN2+μζ7,
with transversality conditions(50)ζjT=0,j∈1,2,3,4,5,6,7,with control functions (n1∗,n2∗) which satisfy the optimality condition given by
(51)n1′t=min1,max0,ζ2−ζ1βI+ηASh1N,n2′t=min1,max0,ζ1−ζ7Sh2.Proof.
With reference to the Hamiltonian (48), the adjoint system (49) is determined by partially differentiating the Hamiltonian (48) with respect to the corresponding state variables S,E,Q,I,A,H,R as
(52)dζ1dt=−∂Hf∂S,dζ2dt=−∂Hf∂E,dζ3dt=−∂Hf∂Q,dζ4dt=−∂Hf∂I,dζ5dt=−∂Hf∂A,dζ6dt=−∂Hf∂H,dζ7dt=−∂Hf∂R.
The characterization of the controls of (51) are derived by solving n1∗ and n2∗ from the equation below.
(53)∂Hf∂n1=0,∂Hf∂n2=0.
Applying bounds on the controls by standard argument, we deduce the characterization.(54)ni∗=0ifωi∗≤0,wi∗if0≤ωi∗≤1,1ifωi∗≥1,where
(55)ω1∗=ζ2−ζ1βI+ηASh1N,ω2∗=ζ1−ζ7Sh2.
Hence. the proof is complete.


## 4.1. Existence of Optimal Controls

As noted in [72], the existence result of Fleming and Richel is considered to show the existence of optimal control duple that minimizes (40) subject to the system (38).Theorem 6.
An optimal control duplen∗ exists that minimizes the objective functional (40) subject to the system (38), given that the below properties are met.
(a)
The set of control is convex and closed(b)
The system (38) is bounded by a linear function in both the state and control variable(c)
The objective functional (40) integrand is convex with respect to the control(d)
There exist constantsb1, b2≥0, and b3≥1 that make the objective functional (40) integrand bounded by the below quantity(42)b1∑i=12ni2b3/2−b2.Proof.
(a)
It is sufficient to writen=n1×n2 by definition of (39). Hence, n=n1×n2 is bounded and convex ∀t∈0,T. In addition, we choose k,d∈N, so that k=k1,k2 and d=d1,d2. Then, ∀v∈0,1, and we have vki+1−vdi∈N, satisfying the convexity property of the control set(b)
We denote the right side of system (38) and the associated solution by f and θ; then,(43)ψ=−1−n1βS1I1N−S2I2N−1−n1βηS1A1N−S2A2N−μS1−S2−n2S1−S21−n1βS1I1N−S2I2N+1−n1βηS1A1N−S2A2N−σE1−E2−τE1−E2−μE1−E2τE1−E2−τ1+μQ1−Q2ασE1−E2−ϕI1−I2−γII1−I2−μI1−I21−ασE1−E2−γAA1−A2−μA1−A2ϕI1−I2+τ1Q1−Q2−δH1−H2−γHH1−H2−μH1−H2γII1−I2+γAA1−A2+n2S1−S2−μR1−R2,where ψ=ft,θ1,n−ft,θ2,n(44)≤2β1−n1S1I1N−S2I2N+2βη1−n1S1A1N−S2A2N+2τE1−E2+2σE1−E2+2ασE1−E2+2τ1Q1−Q2+2γAA1−A2+2γHH1−H2+2n2S1−S2+2ϕI1−I2+μS2−S1+μE2−E1+μQ2−Q1+μI2−I1+μA2−A1+μH2−H1+μR2−R1,≤2β1−n1I1S1−S2+S2I1−I2+2βη1−n1A1S1−S2+S2A1−A2+2τ+2σ+2α+2ασ+μE1−E2+2τ1+μQ1−Q2+2γA+μA1−A2+2γH+μH1−H2+2n2+μS1−S2+2ϕ+μI1−I2+μR2−R1,≤S1−S22β1−n1I1+2βη1−n1A1+2β1−n1S2I1−I2+2βη1−n1S2A1−A2+2τ+2σ+2α+2ασ+μE1−E2+2τ1+μQ1−Q2+2γA+μA1−A2+2γH+μH1−H2+2n2+μS1−S2+2ϕ+μI1−I2+μR2−R1,≤2n2+μS1−S2+2τ+2σ+2α+2ασ+μE1−E2+2τ1+μQ1−Q2+2βΛμ1−n1+2ϕ+μI1−I2+2βηΛμ1−n1+2γA+μA1−A2+2γH+μH1−H2+μR2−R1,≤G1S1−S2+G2E1−E2+G3Q1−Q2+G4I1−I2+G5A1−A2+G6H1−H2+G7R1−R2,
whereG1=2n2+μ, G2=2τ+2σ+2ασ+μ, G3=2τ1+μ, G4=2βΛ/μ1−n1+2ϕ+μ, G5=2βηΛ/μ1−n1+2γA+μ, G6=2γA+μ, G7=μ, and G=maxG1,G2,G3,G4,G5,G6,G7.
Hence,f is uniformly Lipschitz continuous.
(c)
The Langrangian defined asLt,z,n is the integrand of the objective functional (40). Thus, we rewrite Lt,z,n in the form(45)Lt,z,n=r1t,z+r2t,z,with
r1t,z=ν1E+ν2I+ν3A and r2t,z=1/2∑j=12hjnj. The convexity of r2t,z, which is a linear combination of the control function 1/2∑j=12hjnj, needs to be proved. We prove the convexity by letting c:0,12⟶R be hn=1/2n2. Then, ∀u1,u2∈0,12 and ρ∈0,1. Hence, it follows that the below inequality holds.
(46)ρcu1+1−ρcu2≥cρu1+1−ρu2.
This confirms the convexity of the Langrangian with respect to the control.(d)
When observed from (45), we easily see that Lt,z,n≥r1t,z. Hence, we conclude that(47)Lt,z,n≥12h1n12+12h2n22≥b1∑i=12ni2b3/2−b2,with b1=1/2minh1,h2, b2>0, and b3=2. This completes the proof.


## 4.2. Characterization of Optimal Controls

Pontryagin’s maximum principle has been the wheel on which the necessary condition for the COVID (38) duple control needs to meet. With the principle, we convert the COVID (38) and the objective functional (40) into a problem of minimizing the Hamiltonian Hf with respect to the controls njt,j=1,2. Hence, the Hamiltonian Hf is given by
(48)Hf=ν1E+ν2I+ν3A+12h1n12+12h2n22+ζ1Λ−1−n1βSIN−1−n1βSηAN−μS−n2S+ζ21−n1βSIN+1−n1βSηAN−σE−τE−μE+ζ3τE−τ1+μQ+ζ4ασE−ϕI−γII−μI+ζ51−ασE−γAA−μA+ζ6ϕI+τ1Q−δH−γHH−μH+ζ7γII+γAA+γHH+n2S−μR.Theorem 7.
With the optimal control duplen1∗,n2∗ satisfying the condition (41), there exist adjoint variables ζi satisfying the adjoint system below.
(49)dζ1dt=ζ1−ζ21−n1βIS+E+I+A+R+SN2+ζ1−ζ21−n1βηAS+E+I+A+R+SN2+ζ1−ζ7n2+μζ1,dζ2dt=−ν1+ζ2−ζ1βSIN2+ζ2−ζ11−n1βηSAN2+ζ2−ζ5σ+ζ2−ζ3τ+ζ5−ζ4ασ+μζ2,dζ3dt=ζ3−ζ6τ1+μζ3,dζ4dt=−ν2+ζ1−ζ21−n1βSS+E+I+A+R+IN2+ζ2−ζ11−n1βηSAN2+ζ4−ζ6ϕ+ζ4−ζ7γI+μζ4,dζ5dt=−ν3+ζ1−ζ21−n1βηAS+E+I+A+R+AN2+ζ2−ζ1βSIN2+ζ5−ζ7γA+μζ5,dζ6dt=ζ6−ζ7γH+δζ6+μζ6,dζ7dt=ζ2−ζ1βSIN2+ζ2−ζ1βηSAN2+μζ7,
with transversality conditions(50)ζjT=0,j∈1,2,3,4,5,6,7,with control functions (n1∗,n2∗) which satisfy the optimality condition given by
(51)n1′t=min1,max0,ζ2−ζ1βI+ηASh1N,n2′t=min1,max0,ζ1−ζ7Sh2.Proof.
With reference to the Hamiltonian (48), the adjoint system (49) is determined by partially differentiating the Hamiltonian (48) with respect to the corresponding state variables S,E,Q,I,A,H,R as
(52)dζ1dt=−∂Hf∂S,dζ2dt=−∂Hf∂E,dζ3dt=−∂Hf∂Q,dζ4dt=−∂Hf∂I,dζ5dt=−∂Hf∂A,dζ6dt=−∂Hf∂H,dζ7dt=−∂Hf∂R.
The characterization of the controls of (51) are derived by solving n1∗ and n2∗ from the equation below.
(53)∂Hf∂n1=0,∂Hf∂n2=0.
Applying bounds on the controls by standard argument, we deduce the characterization.(54)ni∗=0ifωi∗≤0,wi∗if0≤ωi∗≤1,1ifωi∗≥1,where
(55)ω1∗=ζ2−ζ1βI+ηASh1N,ω2∗=ζ1−ζ7Sh2.
Hence. the proof is complete.


## 5. Model Application with Numerical Examples

As credited to Lenhart and Workman [73], the method of forward-backward sweep has been explored extensively by many researchers as in [74–76] to solve the optimality system of optimal control models numerically. Hence, we consider the method to solve the COVID-19 (1) and the control system (38). We design numerical scheme that uses Runge-Kutta’s fourth-order method [70, 77, 78] to solve the model’s optimality system. The optimality system results and the control problem’s state system are numerically obtained by the Runge-Kutta method of order four, with Matlab. The constants ν1=10, ν2=8, ν3=5, h1=5, and h2=10 were used to balance the terms of the objective functional’s equation, and we utilized the parameter values presented in Table 1 to generate the plots of symptomatic infectious, exposed, and asymptomatic infectious. The weight h2>h1, since we assume that the cost of vaccinating the population would be greater than the self-protection strategy.


### 5.1. Strategy A (withn1 and n2)

In an effort to curtail the pandemic, the controlsn1 and n2 were utilized. The plots of the graphs of Figures 1(a)–1(c) showed an increase in the number of exposed, symptomatic, and asymptomatic infectious individuals in the first 20, 30, and 30 days, respectively, for the noncontrol case. The noncontrol exposed curve rises sharply in the first 5 days until about the 20 days, where its dynamic changes gradually decrease. In a like manner, the symptomatic and asymptomatic infectious curves quickly rise in the early days until 40 days when their dynamics change and begin to decrease. Notwithstanding, utilizing the controls n1 and n2, we notice that the number of exposed, symptomatic, and asymptomatic infectious individuals is greatly minimized. Furthermore, the exposed, symptomatic graphs drop gradually for the first 20 until they completely wipe out the population in 80 and 120 days, respectively. This shows the optimal control strategy’s effectiveness as it has a substantial effect of drastically minimizing the exposed, symptomatic, and asymptomatic graphs. Figure 1(d) depicts the strategy A’s control profile. We noticed that the personal protection control n1 and vaccination control n2 stayed at the upper bound throughout the simulated time of 180 days.Figure 1
Numerical solutions.
(a)(b)(c)(d)


### 5.2. Strategy B (withn1)

Owing to the effort to curtail the disease from spreading, we utilized controlsn1. The plot of the graphs of Figures 2(a)–2(c) of the noncontrol model indicated an increase in the number of exposed, symptomatic, and asymptomatic infectious individuals in the first 10, 35, and 30 days of the simulated time. The exposed graph of (a) quickly rises in the first 10 days until it reaches a maximum height of about 800 when it suddenly changes its dynamics and begins to decrease. The symptomatic and asymptomatic infectious graphs similarly rise in the early days of the simulated time and drop after 35 days. However, with the optimal control strategy of n1, the desired result of minimizing the exposed, symptomatic, and asymptomatic infectious individuals is obtained. Thus, in a similar pattern, even though the control graphs rose early and reached the height of the noncontrol plots, they were lowered than the noncontrol graphs in the final run. This means that the strategy of the controls n1 is efficient for preventing the further spread of the disease, even though it does not bring out the best result. In Figure 2(d), we have a clear view of strategy B’s control profile. As noticed, the control profile of the personal protection control n1 stayed at the middle throughout the simulated time of 180 days.Figure 2
Numerical solutions.
(a)(b)(c)(d)


### 5.3. Strategy C (withn2)

We considered the controln2 in an attempt to defeat the pandemic. The simulated plot of the noncontrol graphs of Figures 3(a)–3(c) indicated an early increase in the number of the exposed, symptomatic, and asymptomatic individuals at the estimated time of about 20, 30, and 30 days, respectively. The situation is reversed with the application of an optimal control strategy. With n2, the plots produced results of the exposed, symptomatic, and asymptomatic individuals’ graphs greatly minimized. Even though the asymptomatic individual’s optimal control plot rises quickly, similar to the noncontrol plot, the control strategy is considered efficient as it minimizes the asymptomatic infectious individuals substantially. Figure 3(d) shows the strategy C’s control profile. We observed that the vaccination control n1 was at the upper bound throughout the simulated time of 180 days.Figure 3
Numerical solutions.
(a)(b)(c)(d)


### 5.4. Conclusion

This research article presented a SEQIAHR compartmental model of COVID-19 to provide insight into the disease’s dynamics by utilizing tailored strategies to minimize the pandemic. We first studied the COVID nonlinear model’s dynamic behaviour by calculating the reproduction number and examining the two nonnegative equilibria’s existence. Global stability analyses for the two equilibria were also carried out by employing the Castillo-Chavez method and Lyapunov function to investigate the global stability of the disease at the disease-free and endemic equilibrium. We carried out a sensitivity analysis on the model to determine the parameters that have relative effects on theR0. We examined the model system (1) to determine the type of bifurcation that it exhibits. Then, we formulated a new SEQIAHR compartmental optimal control model with time-dependent controls: personal protection and vaccination of susceptible individuals and solved it utilizing Pontryagin’s maximum principle after studying the dynamical behaviour of the noncontrol model. We solved the model numerically by considering different simulation controls’ pairing and examined their effectiveness. The results showed that each optimal control strategy chosen has an incomparable impact on the number of the exposed, symptomatic, and asymptomatic individuals compared to the noncontrol model since they substantially minimize exposed, symptomatic, and asymptomatic infectious individuals. Thus, strategy A considered both personal protection and vaccination control. We noticed that the combined effect of the strategies had a significant impact on the disease by emphatically minimizing the exposed, symptomatic, and asymptomatic infectious individuals. Strategy B considered only the personal protection control in its intervention program. Even though, to some extent, the strategy minimized the exposed, symptomatic, and asymptomatic infectious individuals in the long run. We observed that in the early days of the graphs, the strategy struggled to minimize these individuals. In addition, the exposed, symptomatic, and asymptomatic infectious individuals are not greatly minimized. The results of strategy B mean that resorting to this intervention strategy will not bring out the desired results as individuals may refuse to adhere to the personal protection protocols as directed by stakeholders and may put the entire population at risk of the pandemic. Strategy C employed the vaccination control as the only control for its intervention program. However, we obtained a great result as the exposed, symptomatic, and asymptomatic infectious individuals are substantially minimized. The result showed that using only vaccination as the control intervention could have the same mitigating effect on the disease as employing both personal protection and vaccination strategies with a minimized.


## 5.1. Strategy A (withn1 and n2)

In an effort to curtail the pandemic, the controlsn1 and n2 were utilized. The plots of the graphs of Figures 1(a)–1(c) showed an increase in the number of exposed, symptomatic, and asymptomatic infectious individuals in the first 20, 30, and 30 days, respectively, for the noncontrol case. The noncontrol exposed curve rises sharply in the first 5 days until about the 20 days, where its dynamic changes gradually decrease. In a like manner, the symptomatic and asymptomatic infectious curves quickly rise in the early days until 40 days when their dynamics change and begin to decrease. Notwithstanding, utilizing the controls n1 and n2, we notice that the number of exposed, symptomatic, and asymptomatic infectious individuals is greatly minimized. Furthermore, the exposed, symptomatic graphs drop gradually for the first 20 until they completely wipe out the population in 80 and 120 days, respectively. This shows the optimal control strategy’s effectiveness as it has a substantial effect of drastically minimizing the exposed, symptomatic, and asymptomatic graphs. Figure 1(d) depicts the strategy A’s control profile. We noticed that the personal protection control n1 and vaccination control n2 stayed at the upper bound throughout the simulated time of 180 days.Figure 1
Numerical solutions.
(a)(b)(c)(d)


## 5.2. Strategy B (withn1)

Owing to the effort to curtail the disease from spreading, we utilized controlsn1. The plot of the graphs of Figures 2(a)–2(c) of the noncontrol model indicated an increase in the number of exposed, symptomatic, and asymptomatic infectious individuals in the first 10, 35, and 30 days of the simulated time. The exposed graph of (a) quickly rises in the first 10 days until it reaches a maximum height of about 800 when it suddenly changes its dynamics and begins to decrease. The symptomatic and asymptomatic infectious graphs similarly rise in the early days of the simulated time and drop after 35 days. However, with the optimal control strategy of n1, the desired result of minimizing the exposed, symptomatic, and asymptomatic infectious individuals is obtained. Thus, in a similar pattern, even though the control graphs rose early and reached the height of the noncontrol plots, they were lowered than the noncontrol graphs in the final run. This means that the strategy of the controls n1 is efficient for preventing the further spread of the disease, even though it does not bring out the best result. In Figure 2(d), we have a clear view of strategy B’s control profile. As noticed, the control profile of the personal protection control n1 stayed at the middle throughout the simulated time of 180 days.Figure 2
Numerical solutions.
(a)(b)(c)(d)


## 5.3. Strategy C (withn2)

We considered the controln2 in an attempt to defeat the pandemic. The simulated plot of the noncontrol graphs of Figures 3(a)–3(c) indicated an early increase in the number of the exposed, symptomatic, and asymptomatic individuals at the estimated time of about 20, 30, and 30 days, respectively. The situation is reversed with the application of an optimal control strategy. With n2, the plots produced results of the exposed, symptomatic, and asymptomatic individuals’ graphs greatly minimized. Even though the asymptomatic individual’s optimal control plot rises quickly, similar to the noncontrol plot, the control strategy is considered efficient as it minimizes the asymptomatic infectious individuals substantially. Figure 3(d) shows the strategy C’s control profile. We observed that the vaccination control n1 was at the upper bound throughout the simulated time of 180 days.Figure 3
Numerical solutions.
(a)(b)(c)(d)


## 5.4. Conclusion

This research article presented a SEQIAHR compartmental model of COVID-19 to provide insight into the disease’s dynamics by utilizing tailored strategies to minimize the pandemic. We first studied the COVID nonlinear model’s dynamic behaviour by calculating the reproduction number and examining the two nonnegative equilibria’s existence. Global stability analyses for the two equilibria were also carried out by employing the Castillo-Chavez method and Lyapunov function to investigate the global stability of the disease at the disease-free and endemic equilibrium. We carried out a sensitivity analysis on the model to determine the parameters that have relative effects on theR0. We examined the model system (1) to determine the type of bifurcation that it exhibits. Then, we formulated a new SEQIAHR compartmental optimal control model with time-dependent controls: personal protection and vaccination of susceptible individuals and solved it utilizing Pontryagin’s maximum principle after studying the dynamical behaviour of the noncontrol model. We solved the model numerically by considering different simulation controls’ pairing and examined their effectiveness. The results showed that each optimal control strategy chosen has an incomparable impact on the number of the exposed, symptomatic, and asymptomatic individuals compared to the noncontrol model since they substantially minimize exposed, symptomatic, and asymptomatic infectious individuals. Thus, strategy A considered both personal protection and vaccination control. We noticed that the combined effect of the strategies had a significant impact on the disease by emphatically minimizing the exposed, symptomatic, and asymptomatic infectious individuals. Strategy B considered only the personal protection control in its intervention program. Even though, to some extent, the strategy minimized the exposed, symptomatic, and asymptomatic infectious individuals in the long run. We observed that in the early days of the graphs, the strategy struggled to minimize these individuals. In addition, the exposed, symptomatic, and asymptomatic infectious individuals are not greatly minimized. The results of strategy B mean that resorting to this intervention strategy will not bring out the desired results as individuals may refuse to adhere to the personal protection protocols as directed by stakeholders and may put the entire population at risk of the pandemic. Strategy C employed the vaccination control as the only control for its intervention program. However, we obtained a great result as the exposed, symptomatic, and asymptomatic infectious individuals are substantially minimized. The result showed that using only vaccination as the control intervention could have the same mitigating effect on the disease as employing both personal protection and vaccination strategies with a minimized.



---

*Source: 9491847-2022-01-27.xml*