Hamadjam Abboubakar
email: abboubakarhamadjam@yahoo.fr

Jean C Kamgang
email: jckamgang@gmail.com

Daniel Tieudjo
email: tieudjo@yahoo.com

Ensai

Bifurcation Thresholds and Optimal Control in Transmission Dynamics of Arboviral Diseases

Keywords: Arboviral diseases, Bifurcation, Sensitivity analysis, Optimal control, Pontryagin's Maximum Principle, Efficiency analysis, Cost effectiveness analysis AMS Subject Classification (2010): 37N25, 49J15, 92D30

In this paper, we derive and analyse a model for the control of arboviral diseases which takes into account an imperfect vaccine combined with some other mechanisms of control already studied in the literature. We begin by analyse the basic model without controls. We prove the existence of two disease-free equilibrium points and the possible existence of up to two endemic equilibrium points (where the disease persists in the population). We show the existence of a transcritical bifurcation and a possible saddle-node bifurcation and explicitly derive threshold conditions for both, including defining the basic reproduction number, R 0 , which determines whether the disease can persist in the population or not. The epidemiological consequence of saddle-node bifurcation (backward bifurcation) is that the classical requirement of having the reproduction number less than unity, while necessary, is no longer sufficient for disease elimination from the population. It is further shown that in the absence of disease-induced death, the model does not exhibit this phenomenon. We perform the sensitivity analysis to determine the model robustness to parameter values. That is to help us to know the parameters that are most influential in determining disease dynamics. The model is extended by reformulating the model as an optimal control problem, with the use of five time dependent controls, to assess the impact of vaccination combined with treatment, individual protection and vector control strategies (killing adult vectors, reduction of eggs and larvae). By using optimal control theory, we establish optimal conditions under which the disease can be eradicated and we examine the impact of a possible combined control tools on the disease transmission. The Pontryagin's maximum principle is used to characterize the optimal control. Numerical simulations, efficiency analysis and cost effectiveness analysis show that, vaccination combined with other control mechanisms, would reduce the spread of the disease appreciably, and this at low cost.

Introduction

Arboviral diseases are affections transmitted by hematophagous arthropods. There are currently 534 viruses registered in the International Catalog of Arboviruses and 25% of them have caused documented illness in human populations [START_REF] Chippaux | Généralités sur arbovirus et arboviroses-overview of arbovirus and arbovirosis[END_REF][START_REF] Karabatsos | International Catalogue of Arboviruses, including certain other viruses of vertebrates[END_REF][START_REF] Gubler | Human arbovirus infections worldwide[END_REF]. Examples of those kinds of diseases are Dengue, Yellow fever, Saint Louis fever, Encephalitis, West Nile fever and Chikungunya. A wide range of arboviral diseases are transmitted by mosquito bites and constitute a public health emergency of international concern. For example, Dengue, caused by any of four closely-related virus serotypes (DEN-1-4) of the genus Flavivirus, causes 50-100 million infections worldwide every year, and the majority of patients worldwide are children aged 9 to 16 years [START_REF] Sanofi Pasteur | Dengue vaccine, a priority for global health[END_REF][START_REF]Dengue and severe dengue[END_REF][START_REF]Dengue and dengue haemorhagic fever[END_REF].

The dynamics of arboviral diseases like Dengue or Chikungunya are influenced by many factors such as human and mosquito behaviours. The virus itself (multiple serotypes of dengue virus [START_REF]Dengue and severe dengue[END_REF][START_REF]Dengue and dengue haemorhagic fever[END_REF], and multiple strains of chikungunya virus [START_REF] Moulay | Modélisation et analyse mathématique de systèmes dynamiques en épidémiologie[END_REF][START_REF] Parola | Novel chikungunya virus variant in travelers returning from indian ocean islands[END_REF]), as well as the environment, affects directly or indirectly all the present mechanisms of control [START_REF] Brasseur | Analyse des pratiques actuelles destinées á limiter la propagation d'Aedes albopictus dans la zone sud est de la France et propositions d'amélioration[END_REF][START_REF] Carvalho | Mathematical modeling of dengue epidemic: Control methods and vaccination strategies[END_REF]. Indeed, in the absence of conditions which favour the development of their larvae, eggs of certain Aedes mosquitoes (Aedes albopictus, for example) enter in diapause phenomenon, allowing the eggs to hatch even after two years [START_REF] Paupy | Aedes albopictus, an arbovirus vector: from the darkness to the light[END_REF][START_REF] Sota | Survival time and resistance to desiccation of diapause and nondiapause eggs of temperate aedes (stegomyia) mosquitoe[END_REF]. Taking the case of Aedes mosquitoes for example, the main control method used in many countries continues to be space spraying of insecticide for adult mosquito control. This strategy must be repeated constantly, its cost is high, and its effectiveness is limited. Also, Ae. aegypti, for example, prefers to rest inside houses, so truck or aerial insecticide spraying simply does not reach mosquitoes resting in hidden places such as cupboards [START_REF] Parks | Planning social mobilization and communication for dengue fever prevention and control[END_REF]. The different types of control mechanisms put in place to reduce the proliferation of vectors responsible for the transmission of pathogens such as arboviruses are listed below.. (i) Biological control or "biocontrol" is the use of natural enemies to manage vector populations: introduction of parasites, pathogens and predators to target vectors. for example, effective biocontrol agents include predatory fish that feed on mosquito larvae such as mosquitofish (Gambusia affinis) and some cyprinids (carps and minnows) and killifish. Tilapia also consume mosquito larvae [START_REF] Alcaraz | Life history variation of invasive mosquitofish (gambusia holbrooki) along a salinity gradient[END_REF]. As biological control does not cause chemical pollution, it is considered as a better method for mosquito control by many people. However, there are limitations on employing biological agents for mosquito control. The agent introduced usually has to be substantial in number for giving desirable effect.

(ii) Mechanical control consist at the environmental sanitation measures to reduce mosquito breeding sites, such as the physical management of water containers (e.g. mosquito-proof covers for water storage containers, polystyrene beads in water tanks), better designed and reliable water supplies, and recycling of solid waste such as discarded tyres, bottles, and cans [START_REF] Parks | Planning social mobilization and communication for dengue fever prevention and control[END_REF][START_REF] Dumont | Vector control for the chikungunya disease[END_REF].

(iii) Chemical methods [START_REF] Parks | Planning social mobilization and communication for dengue fever prevention and control[END_REF][START_REF] Dumont | Vector control for the chikungunya disease[END_REF]:

-chemical methods against the mosquito's aquatic stages for use in water containers (larviciding -killing of larvae),

-chemical methods directed against adult mosquitoes, such as insecticide space sprays or residual applications (adulticiding -killing of adult mosquitoes),.

(iv) Personal protection consist at the use of repellents, vaporizers, mosquito coils, and insecticide treated screens, curtains, and bednets (for daytime use against Aedes) [START_REF] Parks | Planning social mobilization and communication for dengue fever prevention and control[END_REF].

The main problem encountered in the implementation of some of these control mechanisms is the preservation of the ecological systems. For example, in the "biocontrol" mechanism, direct introduction of tilapia and mosquitofish into ecosystems around the world have had disastrous consequences [START_REF] Alcaraz | Life history variation of invasive mosquitofish (gambusia holbrooki) along a salinity gradient[END_REF]. Also, the chemical methods can not be applied in continuous times. Some chemical product like Deltamethrin seems to be effective only during a couple of hours [START_REF] Dumont | Vector control for the chikungunya disease[END_REF][START_REF] Bosc | Premier Bilan sur les Impacts des Traitements Anti-moustiques, dans le cadre de la lutte Contre le Chikungunya, sur les Espèces et les Milieux de l'île de la Réunion[END_REF][START_REF] Sofia | Optimal Control and Numerical Optimization Applied to Epidemiological Models[END_REF]. So its use over a long period and continuously, leads to strong resistance of the wild populations of Aedes aegypti [START_REF] Darriet | Insecticides larvicides et adulticides alternatifs pour les opérations de démoustication en france, synthèse bibliographique[END_REF], for example.

For all the diseases mentioned above, only yellow fever has a licensed vaccine. Nevertheless, considerable efforts are made to obtain vaccines for other diseases. In the case of dengue, for example, tests carried out in Asia and Latin America, have shown that the future dengue vaccine will have a efficacy between 30.2% and 77.7%, and this, depending on the serotype [START_REF] Sabchareon | Protective efficacy of the recombinant, live-attenuated, cyd tetravalent dengue vaccine in thai schoolchildren: a randomised, controlled phase 2b trial[END_REF][START_REF] Villar | Efficacy of a tetravalent dengue vaccine in children in latin america[END_REF]. Also, the future dengue vaccine will have an overall efficacy of 60.8% against all forms of the disease in children and adolescents aged 9-16 years who received three doses of the vaccine [START_REF]The new england journal of medicine publie les résultats de l'étude clinique d'efficacité de phase iii du candidat vaccin dengue de sanofi pasteur[END_REF].

As the future vaccines (e.g., dengue vaccine) will be imperfect, it is therefore necessary to combine such vaccines with some control mechanisms cited above, to find the best sufficient combination (in terms of efficacy and costs), which permit to decrease the expansion of these kind of diseases in human communities.

A number of studies have been conducted to study host-vector models for arboviral diseases transmission (see [START_REF] Dumont | Vector control for the chikungunya disease[END_REF][START_REF] Abboubakar | Modeling the dynamics of arboviral diseases with vaccination perspective[END_REF][START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF][START_REF] Aldila | An optimal control problem arising from a dengue disease transmission model[END_REF][START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF][START_REF] Blayneha | Backward bifurcation and optimal control in transmission dynamics of west nile virus[END_REF][START_REF] Cannon | An epidemiology model suggested by yellow fever[END_REF][START_REF] Coutinho | Threshold conditions for a non-autonomous epidemic system describing the population dynamics of dengue[END_REF][START_REF] Cruz-Pacheco | Seasonality and outbreaks in west nile virus infection[END_REF][START_REF] Derouich | Dengue fever: mathematical modelling and computer simulation[END_REF][START_REF] Esteva | Analysis of a dengue disease transmission model[END_REF][START_REF]A model for dengue disease with variable human population[END_REF][START_REF] Feng | Competitive exclusion in a vector-host model for the dengue fever[END_REF][START_REF] Garba | Backward bifurcations in dengue transmission dynamics[END_REF][START_REF] Maidana | Dynamic of west nile virus transmission considering several coexisting avian populations[END_REF][START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF][START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF][START_REF] Poletti | Transmission potential of chikungunya virus and control measures: the case of italy[END_REF][START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF]). Some of these works have been conducted to explore optimal control theory for arboviral disease models (see [START_REF] Aldila | An optimal control problem arising from a dengue disease transmission model[END_REF][START_REF] Blayneha | Backward bifurcation and optimal control in transmission dynamics of west nile virus[END_REF][START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF][START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF][START_REF] Weverton | A multiobjective optimization approach for combating aedes aegypti using chemical and biological alternated step-size control[END_REF]).

In [START_REF] Aldila | An optimal control problem arising from a dengue disease transmission model[END_REF], Dipo Aldila and co-workers derive a optimal control problem for a host-vector Dengue transmission model, in which treatments with mosquito repellent are given to adults and children and those who undergo treatment are classified in treated compartments. The only control considered by the authors is the treatment of people with clinical signs of the disease. Blayneh et al. in [START_REF] Blayneha | Backward bifurcation and optimal control in transmission dynamics of west nile virus[END_REF] consider a deterministic model for the transmission dynamics of West Nile virus (WNV) in the mosquito-bird-human zoonotic cycle. They use two control functions, one for mosquito-reduction strategies and the other for personal (human) protection, and redefining the demographic parameters as density-dependent rates. In [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF], Moulay et al. derive optimal prevention (individual protection), vector control (Larvae reduction) and treatment strategies used during the Chikungunya Réunion Island epidemic in 2006. Authors in [START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF] derive the optimal control efforts for vaccination in order to prevent the spread of a Dengue disease using a system of ordinary differential equations (ODEs) for the host and vector populations. Recently, Dias et al. in [START_REF] Weverton | A multiobjective optimization approach for combating aedes aegypti using chemical and biological alternated step-size control[END_REF] analyse the Dengue vector control problem in a multiobjective optimization approach, in which the intention is to minimize both social and economic costs, using a dynamic mathematical model representing the mosquitoes' population. This multiobjective optimization approach consists in finding optimal alternated step-size control policies combining chemical (via application of insecticides) and biological control (via insertion of sterile males produced by irradiation).

None of the above mentioned models [START_REF] Aldila | An optimal control problem arising from a dengue disease transmission model[END_REF][START_REF] Blayneha | Backward bifurcation and optimal control in transmission dynamics of west nile virus[END_REF][START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF][START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF][START_REF] Weverton | A multiobjective optimization approach for combating aedes aegypti using chemical and biological alternated step-size control[END_REF] takes into account the combination of optimal control mechanisms such as vaccination, individual protection, treatment and vector control strategies. In our effort, we investigate such optimal strategies for vaccination combined with individual protection, treatment and two vector controls (adulticiding-killing of adult vectors, and larviciding-killing eggs and larvae), using two systems of ODEs which consist of a complete stage structured model Eggs-Larvae-Pupae for the vectors, and a SEI/SEIR type model for the vector/host population. This provides a new different mathematical perspective to the subject. Furthermore, a efficiency analysis and cost effectiveness analysis, are performed here in order to evaluate the control combination that is most effective in the design of optimal strategies.

We start with the formulation of a model without control which is an modified of the previous models developed in [START_REF] Abboubakar | Modeling the dynamics of arboviral diseases with vaccination perspective[END_REF][START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF]. We compute the net reproductive number N , as well as the basic reproduction number, R 0 , and investigate the existence and stability of equilibria.

We prove that the trivial equilibrium is globally asymptotically stable whenever N < 1. When N > 1 and R 0 < 1, we prove that the system exhibit the backward bifurcation phenomenon. The implication of this occurrence is that the classical epidemiological requirement for effective eradication of the disease, R 0 < 1, is no longer sufficient, even though necessary. We show the existence of a transcritical bifurcation and a possible saddle-node bifurcation and explicitly derive threshold conditions for both.

Then, we formulate an optimal control model by adding five control functions: three for human (vaccination, protection against mosquitoes bites and treatment), and two for mosquitoreduction strategies (the use of adulticide to kill adult vectors, and the use of larvicide to increase the mortality rate of eggs and larvae). By Using optimal control theory, we derive the conditions under which it is optimal to eradicate the disease and examine the impact of a possible combination of vaccination, treatment, individual protection and vector control strategies on the disease transmission. The Pontryagin's maximum principle is used to characterize the optimal control. Numerical simulations, efficiency analysis, as well as, the cost effectiveness analysis, are performed to determine the best combination (in terms of efficacy and cost).

The rest of the paper is organized as follows. In section 2, we present the basic transmission model and carry out some analysis by determining important thresholds such as the net reproductive number N and the basic reproduction number R 0 , and different equilibria of the model. We then demonstrate the stability of equilibria and carry out bifurcation analysis, by deriving the threshold conditions for saddle-node bifurcation. In Section 3 we present the optimal control problem and its mathematical analysis. Section 4 is devoted to numerical simulations, efficiency analysis and cost effectiveness analysis. A conclusion round up the paper.

The basic model and its analysis

The model we propose here is based on the modelling approach given in [START_REF] Abboubakar | Modeling the dynamics of arboviral diseases with vaccination perspective[END_REF][START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF]. For the reader convenience, we briefly recall here main results which are developed in this work.

It is assumed that the human and vector populations are divided into compartments described by time-dependent state variables. The compartments in which the populations are divided are the following ones:

(i) For humans, we consider a SEIR model: Susceptible (denoted by S h ), exposed (E h ), infectious (I h ) and resistant or immune (R h ) which includes naturally-immune individuals. The recruitment in human population is at the constant rate Λ h , and newly recruited individuals enter the susceptible compartment S h . Are concern by recruitment people that are totally naive from the disease. Each human compartment, individual goes out from the dynamics at natural mortality rates µ h . The human susceptible population is decreased following infection, which can be acquired via effective contact with an exposed or infectious vector at a rate

λ h = aβ hv (η v E v + I v ) N h , ( 1
)
where a is the biting rate per susceptible vector, β hv is the transmission probability from an exposed/infectious vector (E v or I v ) to a susceptible human (S h ). The expression of λ h is obtained as follows. The probability that a vector chooses a particular human or other source of blood to bite can be assumed as 1 N h . Thus, a human receives in average a N v N h bites per unit of times. Then, the infection rate per susceptible human is given

aβ hv N v N h (η v E v + I v ) N v .
In expression of λ h , the modification parameter 0 < η v < 1 accounts for the assumed reduction in transmissibility of exposed mosquitoes relative to infectious mosquitoes [START_REF] Abboubakar | Modeling the dynamics of arboviral diseases with vaccination perspective[END_REF][START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF][START_REF] Garba | Backward bifurcations in dengue transmission dynamics[END_REF] (see the references therein for the specific sources). Latent humans (E h ) become infectious (I h ) at rate γ h . Infectious humans recover at a constant rate, σ or dies as consequence of infection, at a disease-induced death rate δ. Immune humans retain their immunity for life. We denote the total human population by N h ,

N h = S h + E h + I h + R h . (2)
(ii) Following [START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF], the stage structured model is used to describe the vector population dynamics, which consists of three main stages: embryonic (E), larvae (L) and pupae (P). Even if eggs (E) and immature stages (L and P) are both aquatic, it is important to dissociate them because, for optimal control point of view, drying the breeding sites does not kill eggs, but only larvae and pupae. Moreover, chemical interventions on the breeding sites have a more great impact on the larval population, but not on the eggs [START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF]. The number of laid eggs is assumed proportional to the number of females. The system of stage structured model of aquatic phase development of vector is given by (see [START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF] for details)

Ė = µ b 1 - E Γ E (S v + E v + I v ) -(s + µ E )E (3a) L = sE 1 - L Γ L -(l + µ L )L (3b) Ṗ = lL -(θ + µ P )P. (3c)
Unlike the authors of [START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF], we take into account the pupal stage in the development of the vector. This is justified by the fact that they do not feed during this transitional stage of development, as they transform from larvae to adults. So, the control mechanisms can not be applied to them. With a rate θ, pupae become female adults. Each vector compartment, individuals goes out from the dynamics at natural mortality rates µ v . The vector susceptible population is decreased following infection, which can be acquired via effective contact with an exposed or infectious human at a rate

λ v = aβ vh (η h E h + I h ) N h , (4)
where β vh is the transmission probability from an exposed/infectious human (E h or I h ) to a susceptible vector (S v ). As well as in the expression of λ h , the modification parameter 0 < η h < 1 in the expression of λ v accounts for the assumed reduction in transmissibility of exposed humans relative to infectious humans [START_REF] Abboubakar | Modeling the dynamics of arboviral diseases with vaccination perspective[END_REF][START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF][START_REF] Garba | Backward bifurcations in dengue transmission dynamics[END_REF]. Latent vectors (E v ) become infectious (I v ) at rate γ v . The vector population does not have an immune class, since it is assumed that their infectious period ends with their death [START_REF]A model for dengue disease with variable human population[END_REF]. So, we denote the total adult vector population by N v ,

N v = S v + E v + I v . (5)

Hosts

Adult vectors

Immature vectors Therefore, our basic arboviral disease model reads as

S h E h I h R h S v E v I v P L E
Ṡh = Λ h -(λ h + µ h )S h (6a) Ėh = λ h S h -(µ h + γ h )E h (6b) İh = γ h E h -[µ h + δ + σ] I h (6c) Ṙh = σI h -µ h R h (6d) Ṡv = θP -λ v S v -µ v S v (6e) Ėv = λ v S v -(µ v + γ v )E v (6f) İv = γ v E v -(µ v )I v (6g) Ė = µ b 1 - E Γ E N v -(s + µ E )E (6h) L = sE 1 - L Γ L -(l + µ L )L (6i) Ṗ = lL -(θ + µ P )P (6j)
where the upper dot denotes the time derivative and λ h and λ v are given by ( 1) and (4), respectively. A schematic of the model is shown in Figure 1. The states and parameters are all strictly positive constants and are described in Table 1 and 2, respectively.

= η v = 0.
(ii) The model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF] is the same that we studied in a previous work (see model 18 in [START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF]) to show that the backward bifurcation is caused by the disease-induced death in human. In this previous work, we have just showed that the occurrence of the backward bifurcation is possible in the model without vaccination. In the present work, we give a sufficient and necessary condition , as well as the explicit expressions of the thresholds which governing this phenomenon.

Basic properties and equilibria

The rates of change of the total populations of humans (2) and adult vectors [START_REF] Brasseur | Analyse des pratiques actuelles destinées á limiter la propagation d'Aedes albopictus dans la zone sud est de la France et propositions d'amélioration[END_REF] for the basic arboviral model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF] are,

Ṅh = Λ h -µ h N h -δI h , Ṅv = θP -µ v N v .
Therefore, by standard arguments (see [START_REF] Abboubakar | Modeling the dynamics of arboviral diseases with vaccination perspective[END_REF][START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF][START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF]) it follows that the feasible region for model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF] is

D = (S h , E h , I h , R h , S v , E v , I v , E, L, P ) ∈ R 10 + : N h ≤ Λ h /µ h ; E ≤ K E ; L ≤ K L ; P ≤ lK L k 7 ; N v ≤ θlK L k 7 k 8 ,
where R 10 + represents the non-negative orthant of R 10 . The model is epidemiologically (the state variables have a valid physical interpretation) and mathematically (the system of equations has a unique solution that is bounded and exists for all time) well-posed in the region D.

For easier readability, we introduce the following quantities,

k 1 := µ h ; k 3 := µ h + γ h ; k 4 := µ h + δ + σ, k 5 := s + µ E ; k 6 := l + µ L ; k 7 := θ + µ P ; k 8 := µ v ; k 9 := µ v + γ v ; k 10 = η h k 4 + γ h ; k 11 = η v k 8 + γ v , (7)
and (the positive quantity),

k 2 = k 3 k 4 -δγ h = µ h k 4 + γ h (µ h + σ).
Without disease in the both populations (i.e

λ h = λ v = 0 or E h = I h = E v = I v = 0
), the basic arboviral model ( 6) have two disease-free equilibria given by E 0 = (N 0 h , 0, 0, 0, 0, 0, 0, 0, 0, 0) which correspond to the trivial equilibrium, and E 1 = (N 0 h , 0, 0, 0, N 0 v , 0, 0, E, L, P ) which correspond to the biological disease-free equilibrium, where

N 0 h = Λ h µ h , N 0 v = Γ E Γ L k 5 k 6 (N -1) µ b (Γ E s + k 6 Γ L ) , P = Γ E Γ L k 5 k 6 k 8 (N -1) µ b θ (Γ E s + k 6 Γ L ) L = Γ E Γ L k 5 k 6 k 7 k 8 (N -1) µ b θl (Γ E s + k 6 Γ L ) , E = Γ E Γ L k 5 k 6 k 7 k 8 (N -1) s (µ b lΓ L θ + k 5 k 7 k 8 Γ E ) , (8)
and N is the net reproductive number [START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF][START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF] given by

N = µ b θls k 5 k 6 k 7 k 8 . (9)
Define the basic reproductive number [START_REF] Diekmann | Mathematical Epidemiology of Infectious Diseases. Model building, analysis and interpretation[END_REF][START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF]]

R 0 = a 2 β hv β vh (γ h + k 4 η h )(γ v + k 8 η v )N 0 v k 3 k 4 k 8 k 9 N 0 h . (10)
We note that,

R 0 = K vh K hv ,
where,

K vh = (K E h vh + K I h vh ) = (a) (β vh ) N 0 v N 0 h 1 k 3 (η h ) + γ h k 4 = aβ vh (γ h + k 4 η h )N 0 v k 3 k 4 N 0 h
, is the number of vector that one human infects through his/her latent/infectious life time. It is equal to the sum of the number of vector infections generated by an exposed human (near the DFE, E 1 ) K E h vh , and the number of vector infections generated by an infectious human (near the DFE) K I h vh . K E h vh is given by the product of the infection rate of exposed humans (aβ vh η h N 0 v /N 0 h ) and the average duration in the exposed (E h ) class (1/k 3 ). K I h vh is given by the product of the infection rate of infectious humans (aβ vh N 0 v /N 0 h ), the probability that an exposed human survives the exposed stage and move to the infectious stage (γ h /(µ h + γ h )) and the average duration in the infectious stage (1/(µ h + δ + σ)). Analogously, we have

K hv = K Ev hv + K Iv hv = (a) (β hv ) 1 k 9 η v + γ v k 8 = aβ hv (γ v + k 8 η h ) k 8 k 9 ,
which is the number of humans that one vector infects through its infectious life time. It is equal to the sum of the number of human infections generated by an exposed vector (near the DFE, E 1 ),K Ev hv , and the number of human infections generated by an infectious vector (near the DFE), K Iv hv . K Ev hv is given by the product of the infection rate of exposed vectors (aβ vh η h ) and the average duration in the exposed (E v ) class (1/(µ v + γ v )). K Iv vh is given by the product of the infection rate of infectious humans (aβ vh ), the probability that an exposed human survives the exposed stage and move to the infectious stage (γ h /k 9 ) and the average duration in the infectious stage (1/k 8 ). The basic reproduction number is equal to the geometric mean of K vh and K hv because infection from human to human goes through one generation of vectors.

The local asymptotic stability result of equilibria E 0 and E 1 is given in the following.

Theorem 1.

(i) if N ≤ 1, the trivial equilibrium E 0 is locally asymptotically stable in D;

(ii) if N > 1, the trivial equilibrium is unstable and the disease-free equilibrium E 1 is locally asymptotically stable in D whenever R 0 < 1.

Proof. See appendix A. .

The epidemiological implication of item (ii) in Theorem 1 is that, in general, when the basic reproduction number, R 0 is less than unity, a small influx of infectious vectors into the community would not generate large outbreaks, and the disease dies out in time (since the DFE is locally asymptotically stable) [START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF][START_REF] Garba | Backward bifurcations in dengue transmission dynamics[END_REF][START_REF] Diekmann | Mathematical Epidemiology of Infectious Diseases. Model building, analysis and interpretation[END_REF][START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF][START_REF] Cruz-Pacheco | Modelling the dynamics of west nile virus[END_REF]. However, we will show in the subsection 2.2 that the disease may still persist even when R 0 < 1.

The global stability of the trivial equilibrium is given by the following result:

Theorem 2. If N ≤ 1
, then E 0 is globally asymptotically stable on D.

Proof. See appendix B.

We now show the existence of endemic equilibria, that is, steady states of model ( 6) where all state variables are positive. First we introduce:

ψ = k 10 aµ h β vh -δγ h k 8 , (11)
R c = 2k 8 k 2 + k 10 aµ h β vh k 3 k 4 k 8 , (12)
R 1b = 1 k 3 k 4 1 k 8 k 9 δγ h (aµ h β vh k 10 + k 2 k 8 ) -(-k 2 ψ) , (13)
R 2b = 1 k 3 k 4 1 k 8 k 9 δγ h (aµ h β vh k 10 + k 2 k 8 ) + (-k 2 ψ) . (14)
Note that (as shown in the Appendix C), when R c < R 0 < 1, ψ ≤ 0 and correspondingly, R 1b and R 2b are real.

With the inequalities,

R c < R 0 < min(1, R 1b ), (15a) max(R c , R 2b ) < R 0 < 1, (15b)
we claim the following result:

Theorem 3. The number of endemic equilibrium points of the basic arboviral disease model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF] depends on R 0 as follows:

(i) For R 0 > 1, the system has a unique endemic equilibrium point.

(ii) For R 0 = 1, the system has (a) A unique endemic equilibrium point if R c < 1.

(b) No endemic equilibrium points otherwise.

(iii) For R 0 < 1, the system has (a) Two endemic equilibrium points if either inequality (15a) or (15b) is satisfied.

(b) A unique endemic equilibrium point if R c < R 0 and either R 0 = R 1b or R 0 = R 2b .
(c) No endemic equilibrium points otherwise.

Proof. See appendix C.

It is clear that case (iii) (item (a)) of theorem 3 indicate the possibility of backward bifurcation (where the locally-asymptotically stable DFE co-exists with a locally asymptotically stable endemic equilibrium when R 0 < 1) in the model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF]. In a previous work (see model 18 in [START_REF] Abboubakar | Backward bifurcation and control in transmission dynamics of arboviral diseases[END_REF]), we just showed that the model exhibited the backward bifurcation phenomenon. In the following, we provide not only a sufficient condition, but also the thresholds which governing this phenomenon.

Bifurcation analysis

Here, we use the centre manifold theory [START_REF] Guckenheimer | Dynamical Systems and Bifurcations of Vector Fields[END_REF] to explore the possibility of bifurcation in (3) at criticality (i.e. the existence and stability of the equilibrium points bifurcating from E 1 at R 0 = 1) by studying the centre manifold near the criticality through the approach developed in [START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF][START_REF] Dushoff | Backward bifurcations and catastrophe in simple models of fatal diseases[END_REF][START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF][START_REF] Carr | Applications of Centre Manifold Theory[END_REF], which is based on general centre manifold theory [START_REF] Guckenheimer | Dynamical Systems and Bifurcations of Vector Fields[END_REF]. To do so, a bifurcation parameter β * hv is chosen, by solving for β hv from R 0 = 1, giving

β * hv = k 3 k 4 k 8 k 9 N 0 h a 2 β vh k 10 k 11 N 0 v . (16)
Let J β * hv denotes the Jacobian of the system (3) evaluated at the DFE (E 1 ) and with

β hv = β * hv . Thus, J =                     -k 1 0 0 0 0 -β * hv η v -β * hv 0 0 0 0 -k 3 0 0 0 β * hv η v β * hv 0 0 0 0 γ h -k 4 0 0 0 0 0 0 0 0 0 σ -µ h 0 0 0 0 0 0 0 - β vh η h S 0 v N 0 h - β vh S 0 v N 0 h 0 -k 8 0 0 0 0 θ 0 β vh η h S 0 v N 0 h β vh S 0 v N 0 h 0 0 -k 9 0 0 0 0 0 0 0 0 0 γ v -k 8 0 0 0 0 0 0 0 K 1 K 1 K 1 -K 2 0 0 0 0 0 0 0 0 0 K 3 -K 4 0 0 0 0 0 0 0 0 0 l -k 7                     , with K 1 = µ b 1 - E * K E , K 2 = k 5 + µ b K E S 0 v . K 3 = s 1 - L * K L , and
K 4 = k 6 + sE * K L .
Note that the system (6), with β hv = β * hv , has a hyperbolic equilibrium point, i.e., the linearised system (6) has a simple eigenvalue with zero real part and all other eigenvalues have negative real part (this follows from the loss of stability of the disease-free equilibrium, E 1 through the transcritical bifurcation). Hence, the centre manifold theory [START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF][START_REF] Dushoff | Backward bifurcations and catastrophe in simple models of fatal diseases[END_REF][START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF][START_REF] Guckenheimer | Dynamical Systems and Bifurcations of Vector Fields[END_REF][START_REF] Carr | Applications of Centre Manifold Theory[END_REF] can be used to analyse the dynamics of the model (3) near β hv = β * hv . The technique in Castillo-Chavez and Song (2004) [START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF] entails finding the left and right eigenvectors of the linearised system above as follows.

The left eigenvector components of J β * hv , which correspond to the uninfected states are zero (see Lemma 3 in [59]). Thus a non-zero components correspond to the infected states. It follows that the matrix J β * hv has a left eigenvector given by v = (v 1 , v 2 , . . . , v 10 ), where

v 1 = v 4 = v 5 = v 8 = v 9 = v 10 = 0; v 2 = k 8 β * hv v 7 ; v 3 = β vh S 0 v (η v k 8 + γ v ) k 4 k 9 N 0 h v 7 ; v 6 = (η v k 8 + γ v ) k 9 v 7 , v 7 > 0. (17)
Similarly, the component of the right eigenvector w are given by w 7 > 0, w 10 > 0,

w 1 = - β * hv (η v k 8 + γ v ) γ v k 1 w 7 ; w 4 = β * hv γ h σ(η v k 8 + γ v ) µ h γ v k 3 k 4 w 7 ; w 2 = µ h k 4 γ h σ w 4 ; w 3 = µ h σ w 4 ; w 5 = - k 8 + γ v γ v w 7 + k 7 K 2 K 4 lK 1 K 3 w 10 ; w 6 = k 8 γ v w 7 ; w 8 = k 7 K 4 lK 3 w 10 ; w 9 = k 7 l w 10 . (18)
Theorem 4.1 in Castillo-Chavez and Song [START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF] is then applied to establish the existence of backward bifurcation in [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF]. To apply such a theorem, it is convenient to let f k represent the right-hand side of the k th equation of the system (3) and let x k be the state variables whose derivative is given by the k th equation for k = 1, . . . , 10. The local bifurcation analysis near the bifurcation point (β hv = β * hv ) is then determined by the signs of two associated constants, denoted by A 1 and A 2 , defined by

A 1 = 10 k,i,j=1 v k w i w j ∂ 2 f k (0, 0) ∂x i ∂x j and A 2 = 10 k,i=1 v k w i ∂ 2 f k (0, 0) ∂x i ∂φ (19)
with φ = β hvβ * hv . It is important to note that in f k (0, 0), the first zero corresponds to the disease-free equilibrium, E 1 , for the system [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF]. Since β hv = β * hv is the bifurcation parameter, it follows from φ = β hvβ * hv that φ = 0 when β hv = β * hv which is the second component in f k (0, 0). Using Eqs. ( 17) and [START_REF] Cruz-Pacheco | Seasonality and outbreaks in west nile virus infection[END_REF] in Eq. ( 19), we obtain

A 1 = v 2 10 i,j=1 w i w j ∂ 2 f 2 (0, 0) ∂x i ∂x j + v 6 10 i,j=1 w i w j ∂ 2 f 6 (0, 0) ∂x i ∂x j and A 2 = v 2 10 i=1 w i ∂ 2 f 2 (0, 0) ∂x i ∂φ (20)
It follows then, after some algebraic computations (see the details in appendix D), that

A 1 = ζ 1 -ζ 2 , (21)
where we have set (see the appendix D for details on derivation of this quantity)

ζ 1 = 2 k 7 K 2 K 4 lK 1 K 3 aβ vh N 0 h (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 v 6 ζ 2 = 2 aβ * hv N 0 h (η v w 6 + w 7 ) (w 2 + w 3 + w 4 )v 2 + 2 aβ vh N 0 h S 0 v N 0 h η h w 2 2 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 + (k 8 + γ v ) γ v (η h w 2 + w 3 ) w 7 v 6 (22)
According to ( 17) and ( 18), we have ζ 1 > 0 and ζ 2 > 0.

We then have

A 2 = aS 0 v N 0 h (η h w 6 + w 7 ) v 2 .
Note that the coefficient A 2 is always positive. Thus, using Theorem 4.1 in [START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF], the following result is established. The direct consequence of Theorem 4 is the following.

Corollary 1. If A 1 < 0 (i.e., ζ 1 < ζ 2 )
, then the unique endemic equilibrium point of the basic model ( 6) is locally asymptotically stable whenever R 0 > 1.

The backward bifurcation phenomenon is illustrated by numerical simulation of the model with the following set of parameter values (it should be noted that these parameters are chosen for illustrative purpose only, and may not necessarily be realistic epidemiologically): Λ h = 30,

β hv = 0.008, η h = 0.78, η v = 0.99, δ = 1, σ = 0.01428, β vh = 0.5, γ v = 1/14, Γ E = 10 4 , Γ L = Γ E /2.
All other parameters are as in Table 2. In this case the conditions required by Theorem 3, case (iii), are satisfied, as well as The associated bifurcation diagram is depicted in figure 2. This clearly shows the coexistence of two locally-asymptotically stable equilibria when R 0 < 1, confirming that the model ( 6) undergoes the phenomenon of backward bifurcation.

ζ 1 = 1.0772 × 10 -7 > ζ 2 = 1.0250 × 10 -9 (so A 1 = 1.0669 × 10 -7 > 0) in Theorem 4. Note, in particular, that with this set of parameters, R c = 0.0367 < 1, R 0 = 0.4359 < 1 (so that R c < R 0 < 1). It follows: d 2 = -5.6537×10 -8 < 0, d 1 = 1.1504 × 10 -10 > 0 and d 0 = -2.4857 × 10 -14 < 0, so that d 2 1 -4d 2 d 0 = 7.6134 × 10 -21 > 0. The resulting two endemic equilibria E 2 = (S * h , E * h , I * h , R * h , S * v , E * v , I * v , E,
The occurrence of the backward bifurcation can be also seen in Figure 3. Here, R 0 is less than the transcritical bifurcation threshold R 0 = 1 (R 0 = 0.4359 < 1), but the solution of the (corresponding to the range (0.0028, 0.0390) of β hv ). The notation EE and DFE stand for endemic equilibrium and disease-free equilibrium, respectively. Solid lines represent stable equilibria and dash lines stand for unstable equilibria. model ( 6) can approach either the endemic equilibrium point or the disease-free equilibrium point, depending on the initial condition.

We know from Theorem 4 that a backward bifurcation scenario is possible for model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF]. Here, we characterize the critical value in terms of a single parameter, the transmission rate β hv , at which the saddle-node bifurcation occurs, i.e. the threshold for the appearance of two endemic equilibria (see Figure 2).

We follow the approach given in [START_REF] Safan | Vaccination based control of infections in sirs models with reinfection: special reference to pertussis[END_REF]. Introducing the quantities,

β = [aµ h β vh k 10 + 2k 2 k 8 ] N 0 h a 2 β vh k 10 k 11 N 0 v , (23)
and,

β ± = N 0 h k 3 k 4 k 10 k 11 a 2 N 0 v β vh δγ h (aµ h β vh k 10 + k 2 k 8 ) ± (-k 2 ψ) 2 , (24)
we have the following result (see the Appendix E for the proof). Theorem 5. Assume that ψ < 0, where ψ is given by [START_REF] Alcaraz | Life history variation of invasive mosquitofish (gambusia holbrooki) along a salinity gradient[END_REF]. Then the backward bifurcation phenomenon takes place in the basic model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF] if and only if

βhv < β hv < min(β -, β * hv ) or max( βhv , β + ) < β hv < β * hv . ( 25
)
The previous analysis is in line with the observation made by Wangari et al. in [START_REF] Wangari | Backward bifurcation in epidemic models: Technical note[END_REF] concerning the bifurcation thresholds of epidemiological models.

2.3 Non-existence of endemic equilibria for R 0 < 1 and δ = 0

In this case, we have the following result. Theorem 6. (i) The model (6) without disease-induced death (δ = 0) has no endemic equilibrium when R 0,δ=0 ≤ 1, and has a unique endemic equilibrium otherwise.

(ii) The DFE, E 1 , of model (6) without disease-induced death (δ = 0), is globally asymptotically stable (GAS) in D if R 0,δ=0 < 1.

Proof. See appendix F.

Sensitivity analysis

We carried out the sensitivity analysis to determine the model robustness to parameter values. That is to help us to know the parameters that are most influential in determining disease dynamics. Following the approach by Marino et al. [START_REF] Marino | A methodology for performing global uncertainty and sensitivity analysis in systems biology[END_REF] and Wu et al. [START_REF] Wu | Sensitivity analysis of infectious disease models: methods, advances and their application[END_REF], partial rank correlation coefficients (PRCC) between the basic reproduction number R 0 and each parameter are derived from 5,000 runs of the Latin hypercube sampling (LHS) method, which is a stratified Monte Carlo sampling method that divides each parameter's range into N equal intervals and randomly draws one sample from each interval [START_REF] Wu | Sensitivity analysis of infectious disease models: methods, advances and their application[END_REF][START_REF] Stein | Large sample properties of simulations using latin hypercube sampling[END_REF]. The parameters are assumed to be random variables with uniform distributions with their mean value listed in Table 2.

With these 5,000 runs of LHS, the derived distribution of R 0 is given in Figure 4. This sampling shows that the mean of R 0 is 1.9583 and the standard deviation is 1.8439. This implies that for the mean of parameter values given in Table 2, we may be confident that the model predicts an endemic state, since the basic reproduction number is greater that unity. The probability that R 0 > 1 (the disease-free equilibrium is unstable and there is exactly one endemic equilibrium point) is 64.48%. We also evaluate the probabilities that conditions (i), (ii) and (iii) in Theorem 3 are satisfied. Let us set P[X] the probability of X, and the sets of parameter values for which (N > 1) is true by Φ 1 , the sets of parameter values for which

(R c < R 0 < min(1, R 1b )) or (max(R c , R 2b ) < R 0 < 1) by Φ 2 ,

and the sets of parameter values for which

{(R 0 = R 1b ) or (R 0 = R 2b )} by Φ 3 , P [¬Φ 1 ] = P [N ≤ 1] = 0.0020, (26a)
P [Φ 1 ] = P [N > 1] = 0.9980, (26b)
P [Φ 1 and (R 0 < 1) and Φ 2 ] = 0.0026, (26c)

P [Φ 1 and (R 0 < 1) and Φ 3 ] = 0, (26d)
P [Φ 1 and(R 0 < 1) and ¬Φ 2 and ¬Φ 3 ] = 0.3526, (26e)

P [Φ 1 and (R 0 < 1)] = 0.3552, (26f)
P [Φ 1 and (R 0 ≥ 1)] = 0.6448. ( 26g
)
Therefore, the probability that the trivial disease-free equilibrium is locally asymptotically stable is 0.0020 (from (26a)), the probability that the disease free equilibrium point is locally asymptotically stable is 0.3552 (from (26f)), the probability that the disease free equilibrium point is locally asymptotically stable and (i) there are no endemic equilibrium points is 35.26% (ii) there are exactly one endemic equilibrium point is 0 (from (26d)), (iii) there are exactly two endemic equilibrium points is 0.0026 (from (26c)). This implies that for the ranges of parameter values given in Table 4, the disease-free equilibrium point is likely to be locally asymptotically stable and the probability of co-existence of a locally asymptotically stable endemic equilibrium point (occurrence of backward bifurcation phenomenon) is very small and insignificant.

We now use sensitivity analysis to analyse the influence of each parameter on the basic reproductive number. From the previously sampled parameter values, we compute the PRCC between R 0 and each parameter of model ( 3). The parameters with large PRCC values (> 0.5 or < -0.5) statistically have the most influence [START_REF] Wu | Sensitivity analysis of infectious disease models: methods, advances and their application[END_REF]. The results, displayed in Table 3 and Figure 5, show that the basic reproduction number R 0 outcome measures are sensitive to changes in the parameters β vh , β hv , a, θ, µ v , Λ h and l.

The sensitivity results suggest that the control of the epidemic of arboviral diseases pass through a combination of immunization against arbovirus, individual protection against vector bites, treatment of infected human, and vector control mechanisms.

A Model for Optimal Control

There are several possible interventions in order to reduce or limit the proliferation of mosquitoes and the explosion of the number of infected humans and mosquitoes. In addition of controls used in [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF], we add vaccination and the control of adult vectors as control variables to reduce or even eradicate the disease. So we introduce five time dependent controls:

1. The first control 0 ≤ u 1 (t) ≤ 1 denotes the percentage of susceptible individuals that one decides to vaccinate at time t. A parameter ω associated to the control u 1 (t) represents the waning immunity process [START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF].

2. The second control 0 ≤ u 2 (t) ≤ 1 represents efforts made to protect human from mosquito bites. It mainly consists to the use of mosquito nets or wearing appropriate clothes [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF]. Thus we modify the infection term as follows:

λ c h = (1 -α 1 u 2 (t))λ h , , λ c v = (1 -α 1 u 2 (t))λ v ( 27
)
where α 1 measures the effectiveness of the prevention measurements against mosquito bites.

3. The third control 0 ≤ u 3 (t) ≤ 1 represents efforts made for treatment. It mainly consists in isolating infected patients in hospitals, installing an anti-mosquito electric diffuser in the hospital room, or symptomatic treatments [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF]. Thus we modify the recovery rate such that σ c h := σ h + α 2 u 3 . α 2 is the effectiveness of the anti-arboviral diseases drugs with α 2 = 0.3 [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF]. Note that this control also permit to reduce the disease-induced death.

4. The fourth control 0 ≤ u 4 (t) ≤ 1 represents mosquitoes adulticiding effort with killing efficacy c m . Thus the mosquito natural mortality rate becomes µ c v = µ v + c m u 4 (t). 5. The fifth control 0 ≤ u 5 (t) ≤ 1 represents the effect of interventions used for the vector control. It mainly consists in the reduction of breeding sites with chemical application methods, for instance using larvicides like BTI (Bacillus Thuringensis Israelensis) which is a biological larvicide, or by introducing larvivore fish. This control focuses on the reduction of the number of larvae, and thus eggs, of any natural or artificial water-filled container [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF]. Thus the eggs and Larvae natural mortality rate become µ c E = µ E +η 1 u 5 (t) and µ c L = µ L + η 2 u 5 (t) where η 1 , η 2 , represent the chemical eggs and larvae mortality rate, respectively [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF].

Note that 0 ≤ u i ≤ 1, for i = 1, . . . , 5, means that when the control is zero there is no any effort invested (i.e. no control) and when it is one, the maximum control effort is invested.

Therefore, our optimal control model of arboviral diseases reads as

                                         Ṡh = Λ h -[(1 -α 1 u 2 (t))λ h + µ h + u 1 (t)] S h + ωu 1 (t)R h Ėh = (1 -α 1 u 2 (t))λ h S h -(µ h + γ h )E h İh = γ h E h -[µ h + (1 -α 2 u 3 (t))δ + σ + α 2 u 3 (t)] I h Ṙh = (σ + α 2 u 3 (t))I h + u 1 S h -(µ h + ωu 1 )R h Ṡv = θP -(1 -α 1 u 2 (t))λ v S v -(µ v + c m u 4 (t))S v Ėv = (1 -α 1 u 2 (t))λ v S v -(µ v + γ v + c m u 4 (t))E v İv = γ v E v -(µ v + c m u 4 (t))I v Ė = µ b 1 - E Γ E (S v + E v + I v ) -(s + µ E + η 1 u 5 (t))E L = sE 1 - L Γ L -(l + µ L + η 2 u 5 (t))L Ṗ = lL -(θ + µ P )P (28)
with initial conditions given at t = 0. For the non-autonomous system [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF], the rate of change of the total populations of humans and adults vectors is given, respectively, by

Ṅh = Λ h -µ h N h -(1 -α 2 u 3 (t))δI h Ṅv = θP -(µ v + c m u 4 (t))N v (29)
For bounded Lebesgue measurable controls and non-negative initial conditions, non-negative bounded solutions to the state system exist [START_REF] Lukes | Differential equations : classical to controlled[END_REF].

The objective of control is to minimize: the number of symptomatic humans infected with arboviruses (that is, to reduce sub-population I h ), the number of vector (N v ) and the number of eggs and larvae (that is, to reduce sub-population E and L, respectively), while keeping the costs of the control as low as possible.

To achieve this objective we must incorporate the relative costs associated with each policy (control) or combination of policies directed towards controlling the spread of arboviral diseases. We define the objective function as

J(u 1 , u 2 , u 3 , u 4 , u 5 ) = t f 0 D 1 I h (t) + D 2 N v (t) + D 3 E(t) + D 4 L(t) + 5 i=1 B i u 2 i (t) dt (30)
and the control set

∆ = {(u 1 , u 2 , u 3 , u 4 , u 5 )|u i (t) is Lebesgue measurable on [0, t f ], 0 ≤ u i (t) ≤ 1, i = 1, . . . , 5}.
The first fourth terms in the integrand J represent benefit of I h , N v , E and L populations, describing the comparative importance of the terms in the functional. A high value of D 1 for example, means that it is more important to reduce the burden of disease as reduce the costs related to all control strategies [START_REF] Buonomo | A simple analysis of vaccination strategies for rubella[END_REF]. Positive constants B i , i = 1, . . . , 5 are weight for vaccination, individual protection (human), treatment and vector control effort respectively, which regularize the optimal control. In line with the authors of some studies on the optimal control (see [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF][START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF][START_REF] Weverton | A multiobjective optimization approach for combating aedes aegypti using chemical and biological alternated step-size control[END_REF][START_REF] Adams | Dynamic multidrug therapies for hiv: optimal and sti control approaches[END_REF][START_REF] Buonomo | A simple analysis of vaccination strategies for rubella[END_REF][START_REF] Zaman | Stability analysis and optimal vaccination of an sir epidemic model[END_REF][START_REF] Jung | Optimal control of treatments in a two-strain tuberculosis model[END_REF][START_REF] Yusuf | Optimal control of vaccination and treatment for an sir epidemiological model[END_REF]), we choose a linear function for the cost on infection,

D 1 I h , D 2 N v , D 3 E, D 4
L, and quadratic forms for the cost on the controls

B 1 u 2 1 , B 2 u 2 2 , B 3 u 2 3 , B 4 u 2 4
, and B 5 u 2 5 . This choice can be justified by the following arguments: (i) An epidemiological control can be likened to an expenditure of energy, by bringing to the applications of physics in control theory;

(ii) In a certain sense, minimize u i is like minimize u 2 i , because u i > 0, i = 1, . . . , 5. (iii) The quadratic controls give rise to controls as feedback law, which is convenient for calculations.

We solve the problem using optimal control theory.

Theorem 7. Let X = (S h , E h , I h , R h , S v , E v , I v , E, L, P ). The following set Ω = X ∈ R 10 : N h ≤ Λ h µ h ; E ≤ Γ E ; L ≤ Γ L ; P ≤ lΓ L k 7 ; N v ≤ θlΓ L k 7 k 8
is positively invariant under system (28).

Proof. On the one hand, one can easily see that it is possible to get,

                                     Ṡh ≥ -(λ h + µ h ) S h Ėh ≥ -(µ h + γ h )E h İh ≥ -(µ h + δ + σ)I h Ṙh ≥ -µ h R h Ė ≥ -( µ b K E + s + µ E + η 1 )E L ≥ -( s K L + l + µ L + η 2 )L Ṗ ≥ -(θ + µ P + η 3 )P Ṡv ≥ -(λ v + µ v )S v Ėv ≥ -(µ v + γ v )E v İv ≥ -µ v I v (31)
for (S h (0), E h (0), I h (0), R h (0), E(0), A(0), P (0), S v (0), E v (0), I v (0)) ≥ 0. Thus, solutions with initial value in Ω remain nonnegative for all t ≥ 0. On the other hand, we have

                   Ṅh ≤ Λ h -µ h N h Ṅv ≤ θP -µ v N v Ė ≤ µ b 1 - E K E (S v + E v + I v ) -(s + µ E )E L ≤ sE 1 - L K L -(l + µ L )L Ṗ ≤ lL -(θ + µ P )P (32)
The right hand side of the inequalities correspond to the transmission model without control, and it is easy to show that solutions remain in Ω. Then using Gronwall's inequality, we deduce that solutions of (28) are bounded.

Existence of an optimal control

The existence of an optimal control can be obtained by using a result of Fleming and Rishel [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF].

Theorem 8. Consider the control problem with system (28).

There exists

u ⋆ = (u ⋆ 1 , u ⋆ 2 , u ⋆ 3 , u ⋆ 4 , u ⋆ 5 ) such that min (u 1 ,u 2 ,u 3 ,u 4 ,u 5 )∈∆ J(u 1 , u 2 , u 3 , u 4 , u 5 ) = J(u ⋆ 1 , u ⋆ 2 , u ⋆ 3 , u ⋆ 4 , u ⋆ 5 )
Proof. To use an existence result, Theorem III.4.1 from [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF], we must check if the following properties are satisfied:

1-the set of controls and corresponding state variables is non empty; 2-the control set ∆ is convex and closed;

3-the right hand side of the state system is bounded by a linear function in the state and control;

4-the integrand of the objective functional is convex;

5-there exist constants c 1 > 0 , c 2 > 0 , and β > 1 such that the integrand of the objective functional is bounded below by c 1 5

i=1 |u i | 2 β 2 -c 2 .
In order to verify these properties, we use a result from Lukes [START_REF] Lukes | Differential equations : classical to controlled[END_REF] to give the existence of solutions for the state system [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] with bounded coefficients, which gives condition 1. Since by definition, the control set ∆ is bounded , then condition 2 is satisfied. The right hand side of the state system (28) satisfies condition 3 since the state solutions are bounded. The integrand of our objective functional is clearly convex on ∆, which gives condition 4. There are c 1 > 0,

c 2 > 0 and β > 1 satisfying D 1 I h + D 2 N v + D 3 E + D 4 L + 5 i=1 B i u 2 i ≥ c 1 5 i=1 |u i | 2 β 2 -c 2 ,
because the states variables are bounded. Thus condition 5 is satisfied. We conclude that there exists an optimal control

u * = (u ⋆ 1 , u ⋆ 2 , u ⋆ 3 , u ⋆ 4 , u ⋆ 5
) that minimizes the objective functional J (u 1 , u 2 , u 3 , u 4 , u 5 ).

Characterization of an optimal control

The necessary conditions that an optimal control must satisfy come from the Pontryagin's Maximum Principle (PMP) [START_REF] Pontryagin | The mathematical theory of optimal processes[END_REF]. This principle converts (28)-( 30) into a problem of minimizing point wise a Hamiltonian H, with respect to (u 1 , u 2 , u 3 , u 4 , u 5 ):

H = D 1 I h + D 2 N v + D 3 E + D 4 L + 5 i=1 B i u 2 i + λ S h {Λ h -[(1 -α 1 u 2 )λ h + µ h + u 1 ] S h + ωu 1 R h } + λ E h {[1 -α 1 u 2 ] λ h S h -(µ h + γ h )E h } + λ I h {γ h E h -(µ h + (1 -α 2 u 3 )δ + σ + α 2 u 3 )I h } + λ R h {(σ + α 2 u 3 )I h + u 1 S h -(µ h + ωu 1 )R h } + λ Sv {θP -[1 -α 1 u 2 ] λ v S v -(µ v + c m u 4 )S v } + λ Ev {(1 -α 1 u 2 ) λ v S v -(µ v + γ v + c m u 4 )E v } + λ Iv {γ v E v -(µ v + c m u 4 )I v } + λ E µ b 1 - E Γ E (S v + E v + I v ) -(s + µ E + η 1 u 5 )E + λ L sE 1 - L Γ L -(l + µ L + η 2 u 5 )L + λ P {lL -(θ + µ P )P } (33)
where the λ i , i = S h , E h , I h , R h , S v , E v , I v , E, L, P are the adjoint variables or co-state variables. Applying Pontryagin's Maximum Principle [START_REF] Pontryagin | The mathematical theory of optimal processes[END_REF], we obtain the following result.

Theorem 9. Given an optimal control

u ⋆ = (u ⋆ 1 , u ⋆ 2 , u ⋆ 3 , u ⋆ 4 , u ⋆ 5 ) and solutions (S ⋆ h , E ⋆ h , I ⋆ h , R ⋆ h , S ⋆ v , E ⋆ v , I ⋆ v , E ⋆ , A ⋆ , P ⋆ ) of the corresponding state system (28), there exist adjoint variables Π = (λ S h , λ E h , λ I h , λ R h , λ Sv , λ Ev , λ Iv , λ E , λ L , λ P ) satisfying, dλ S h dt = µ h λ S h + u 1 (λ S h -λ R h ) + (1 -α 1 u 2 )λ h 1 - S h N h (λ S h -λ E h ) + (1 -α 1 u 2 ) S v λ v N h (λ Ev -λ Sv ) (34) dλ E h dt = µ h λ E h + γ h (λ E h -λ I h ) + (1 -α 1 u 2 ) S h λ h N h (λ E h -λ S h ) + (1 -α 1 u 2 ) S v N h (aβ vh η h -λ v ) (λ Sv -λ Ev ) ( 35
)
dλ I h dt = -D 1 + [µ h + (1 -α 2 u 3 )δ] λ I h + (σ + α 2 u 3 )(λ I h -λ R h ) + (1 -α 1 u 2 ) S h λ h N h (λ E h -λ S h ) + (1 -α 1 u 2 ) S v N h (aβ vh -λ v ) (λ Sv -λ Ev ) (36) dλ R h dt = µ h λ R h + ωu 1 (λ R h -λ S h ) + (1 -α 1 u 2 ) S h λ h N h (λ E h -λ S h ) + (1 -α 1 u 2 ) S v λ v N h (λ Ev -λ Sv ) ( 37
)
dλ Sv dt = -D 2 + (µ v + c m u 4 )λ Sv + (1 -α 1 u 2 )λ v (λ Sv -λ Ev ) -µ b 1 - E Γ E λ E ( 38
)
dλ Ev dt = -D 2 + (µ v + c m u 4 )λ Ev + γ v (λ Ev -λ Iv ) + aη v β hv (1 -α 1 u 2 )(λ S h -λ E h ) S h N h -µ b 1 - E Γ E λ E (39) dλ Iv dt = -D 2 + (µ v + c m u 4 )λ Iv + aβ hv (1 -α 1 u 2 ) S h N h (λ S h -λ E h ) -µ b 1 - E Γ E λ E ( 40
)
dλ E dt = -D 3 + µ b Γ E N v + s + µ E + η 1 u 5 λ E -s 1 - L Γ L λ L ( 41
)
dλ L dt = -D 4 -lλ P + s Γ L E + µ L + l + η 2 u 5 λ L ( 42
)
dλ P dt = (µ P + θ)λ Pθλ Sv [START_REF] Moulay | The chikungunya disease: Modeling, vector and transmission global dynamics[END_REF] and the transversality conditions

λ * i (t f ) = 0, i = 1, . . . 10. (44)
Furthermore,

u ⋆ 1 = min 1, max 0, (S h -ωR h )(λ S h -λ R h ) 2B 1 , u ⋆ 2 = min 1, max 0, α 1 [λ h S h (λ E h -λ S h ) + λ v S v (λ Ev -λ Sv )] 2B 2 , u ⋆ 3 = min 1, max 0, α 2 [(1 -δ)λ I h -λ R h ] I h 2B 3 , u ⋆ 4 = min 1, max 0, c m [S v λ Sv + E v λ Ev + I v λ Iv ] 2B 4 , u ⋆ 5 = min 1, max 0, η 1 Eλ E + η 2 Lλ L 2B 5 . (45)
Proof. The differential equations governing the adjoint variables are obtained by differentiation of the Hamiltonian function, evaluated at the optimal control. Then the adjoint system can be written as

dλ S h dt = - ∂H ∂S h , dλ E h dt = - ∂H ∂E h , dλ I h dt = - ∂H ∂I h , dλ R h dt = - ∂H ∂R h , dλ Sv dt = - ∂H ∂S v , dλ Ev dt = - ∂H ∂E v , dλ Iv dt = - ∂H ∂I v , dλ E dt = - ∂H ∂E , dλ L dt = - ∂H ∂L , dλ P dt = - ∂H ∂P ,
with zero final time conditions (transversality).

To get the characterization of the optimal control given by (45), we follow [START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF][START_REF] Lenhart | Optimal Control Applied to Biological Models[END_REF] and solve the equations on the interior of the control set,

∂H ∂u i = 0, i = 1, . . . , 5.
Using the bounds on the controls, we obtain the desired characterization. This ends the proof.

Numerical simulations and discussion

The simulations were carried out using the values of Table 4. We use an iterative scheme to solve the optimality system. We first solve the state equations ( 28) with a guess for the controls over the simulated time using fourth order Runge-Kutta scheme. Then, we use the current iterations solutions of the state equation to solve the adjoint equations ( 34)-( 43) by a backward fourth order Runge-Kutta scheme. Finally, we update the controls by using a convex combination of the previous controls and the value from the characterizations [START_REF] Okosun | Optimal control analysis of a malaria disease transmission model that includes treatment and vaccination with waning immunity[END_REF] (see e.g. [START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF][START_REF] Buonomo | A simple analysis of vaccination strategies for rubella[END_REF][START_REF] Zaman | Stability analysis and optimal vaccination of an sir epidemic model[END_REF][START_REF] Lenhart | Optimal Control Applied to Biological Models[END_REF][START_REF] Okosun | Optimal control analysis of a malaria disease transmission model that includes treatment and vaccination with waning immunity[END_REF]). The values chosen for the weights in the objective functional J (see Eq. ( 30)) are given in Table 5. Table 6 gives the initials conditions of state variables. We simulated the system (28) in a period of twenty days (t f = 20).

As the purpose of our study is to seek the best combination linking vaccination to other control mechanism, we will just determine the best strategy among the strategies listed in Table 7. Therefore, we will distinguished five control strategies at follows: Vaccine combined with individual protection, treatment and larvicide Z 3

Vaccine combined with treatment, adulticide and larvicide Z 4

Vaccine combined with individual protection, adulticide and larvicide Z Combination of the five controls (i) Vaccination combined with individual protection, treatment and adulticide With this strategy, only the combination of the control u 1 on vaccination, the control u 2 on individual protection, the control u 3 on treatment and the control u 4 on adulticide, is used to minimise the objective function J (30), while the other control u 5 is set to zero. On figure 7, we observed that the control strategy resulted in a decrease in the number of infected humans (I h ) while an increase is observed in the number of infected humans (I h ) in strategy without control. The use of this combination have also a great impact on the decreasing total vector population (N v ), as well as aquatic vector populations (E and L).

(ii) Vaccination combined with individual protection, treatment and larvicide

With this strategy, only the combination of the control u 1 on vaccination, the control u 2 on individual protection, the control u 3 on treatment and the control u 5 on larvicide, is used to minimise the objective function J (30), while the other control u 4 is set to zero. On figure 8, we observed that the control strategy resulted in a decrease in the number of infected humans [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with individual protection, treatment and adulticide (u 1 = 0, u 2 = 0, u 3 = 0, u 4 = 0 ). (I h ) while an increase is observed in the number of infected humans (I h ) in strategy without control. The use of this combination have no impact on the decreasing total vector population (N v ), as well as aquatic vector populations (E and L). [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with individual protection, treatment and larvicide (u 1 = 0, u 2 = 0, u 3 = 0, u 5 = 0 ).

(iii) Vaccination combined with treatment, adulticide and larvicide

With this strategy, only the combination of the control u 1 on vaccination, the control u 3 on treatment, the control u 4 on adulticide and the control u 5 on larvicide, is used to minimise the objective function J (30), while the other control u 2 is set to zero. On figure 9, we observed that the control strategy resulted in a decrease in the number of infected humans (I h ) while an increase is observed in the number of infected humans (I h ) in strategy without control. The use of this combination have a considerable impact on the decreasing total vector population (N v ), as well as aquatic vector populations (E and L). [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with treatment, adulticide and larvicide (u 1 = 0, u 3 = 0, u 4 = 0, u 5 = 0 ).

(iv) Vaccination combined with individual protection, adulticide and larvicide

With this strategy, only the combination of the control u 1 on vaccination, the control u 2 on individual protection, the control u 4 on adulticide and the control u 5 on larvicide, is used to minimise the objective function J (30), while the other control u 4 are set to zero. On figure 10, we observed that the control strategy resulted in a decrease in the number of infected humans (I h ) while an increase is observed in the number of infected humans (I h ) in strategy without control. The use of this combination have a great impact on the decreasing total vector population (N v ), as well as aquatic vector populations (E and L).

(v) The combination of all the five controls

In this strategy, the combination of all the five controls is applied. On figure 11, we observed that combining all the five controls give a better result in a decrease in the number of infected humans (I h ), as well as, the total number of vector population (N v ), and the aquatic vector populations (E and L).

Figures 7, 9, 10, and 11 have the same shape, which it is difficult to say what is the best control strategy. In the following, we make an efficiency analysis and a cost effectiveness analysis to determine the best strategy in terms of efficiency and cost [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with individual protection, adulticide and larvicide (u 1 = 0, u 2 = 0, u 4 = 0, u 5 = 0 ). [START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using the combination of all the five controls (u i = 0, i = 1, . . . , 5).

Efficiency analysis

In line with Yang and Ferreira [START_REF] Yang | Assessing the effects of vector control on dengue transmission[END_REF], Dumont and Chiroleu [START_REF] Dumont | Vector control for the chikungunya disease[END_REF], and Carvalho et al. [START_REF] Carvalho | Mathematical modeling of dengue epidemic: Control methods and vaccination strategies[END_REF], we compare the effects of different strategies applied on the arboviral diseases, by the introduction of the efficiency index, designed by F. To this aim, we define the variable A as the area comprised between the curve of the symptomatic infectious human (I h ) population size, for instance, and the time axis during the period of time from 0 to t f , as

A = t h 0 I h (t)dt, (46)
which measures the cumulated number of infectious human during the time interval [0, t f ] [START_REF] Carvalho | Mathematical modeling of dengue epidemic: Control methods and vaccination strategies[END_REF][START_REF] Yang | Assessing the effects of vector control on dengue transmission[END_REF]. Hence the efficiency index, F, be can defined by

F = 1 - A c h A (0) h × 100, ( 47
)
where A c h and A 0 h are the cumulated number of infectious human with and without the different controlling mechanisms, respectively. So, It follows that the best strategy will be the one whom efficiency index will be the biggest [START_REF] Carvalho | Mathematical modeling of dengue epidemic: Control methods and vaccination strategies[END_REF][START_REF] Yang | Assessing the effects of vector control on dengue transmission[END_REF].

With the previous simulations, we resume the efficiency index of different strategies in the Tables 8.

Table 8: Table of efficiency index (the case of infected humans).

Strategy

A

I h = F I h = 100× Strategy A I h = F I h = 100× t h 0 I h (t)dt 1 - A c I h A (0) I h t h 0 I h (t)dt 1 - A c I h A (0) I h
No controls 4.1052 × 10

Cost Effectiveness Analysis

The analysis of efficiency allowed us to determine the most efficient strategy, regardless of the cost associated with each control. In what follows, we will among the strategies listed in Table 7, determine which one is the most efficient and which can be implemented at lower cost. To this aim, we follow Okosun and co-workers [START_REF] Okosun | Optimal control analysis of a malaria disease transmission model that includes treatment and vaccination with waning immunity[END_REF], and use cost effectiveness analysis to determine the most cost effective strategy to use the controls of arboviral diseases (Strategies Z 1 -Z).

To this aim, we calculate the Incremental Cost-Effectiveness Ratio (ICER) which is generally described as the additional cost per additional health outcome (see [START_REF] Okosun | Optimal control analysis of a malaria disease transmission model that includes treatment and vaccination with waning immunity[END_REF]).

Based on the model simulation results, we rank the strategies in decreasing number of Total number of infectious individuals (see Table 9). The difference between the total number of infectious individuals without control and the total number of infectious individuals with control was used to determine the "total number of infection averted" used in the tables of cost-effectiveness analysis.

We obtain the ICER of strategies (i) and (j) by applying the formula given by equation [START_REF] Poletti | Transmission potential of chikungunya virus and control measures: the case of italy[END_REF].

ICER(i) = T otal cost(i) T otal inf ection averted(i) ICER(j) =
T otal cost(j) -T otal cost(i) T otal inf ection averted(j) -T otal inf ection averted(i)

(48)
The comparison between ICER(Z 4 ) and ICER(Z 3 ) shows a cost saving of $ 993.6407 for strategy Z 3 over strategy Z 4 . The negative ICER for strategy Z 3 indicates that the strategy Z 4 is "strongly dominated". That is, strategy Z 4 is more costly and less effective than strategy Z 3 . Therefore, strategy Z 4 , the strongly dominated is excluded from the set of alternatives so it does not consume limited resources. We exclude strategy Z 4 and compare strategy Z 3 with Z 2 . From table 9, we have:

The comparison between ICER(Z 3 ) and ICER(Z 2 ) shows a cost saving of $ 8569.4175 for strategy Z 3 over strategy Z 2 . That is, strategy Z 2 is more costly and less effective than strategy Z 3 . Therefore, strategy Z 2 , the strongly dominated is excluded.

We then compare strategy Z 3 with Z 1 . From table 9, we have:

The comparison between ICER(Z 3 ) and ICER(Z 1 ) shows a cost saving of $ 1333.3333 for strategy Z 1 over strategy Z 3 . So the strategy Z 3 is "strongly dominated". That is, strategy Z 3 is more costly and less effective than strategy Z 1 . Therefore, strategy Z 3 , the strongly dominated is excluded.

We then compare strategy Z 1 with Z. From table 9, we have From Table 13, it follows that strategy Z 3 is equivalent in term of efficiency and cost at strategy Z.

With these results, we conclude that the strategies Z (combination of the five control) and Z 3 (vaccination u 1 combined with personal protection u 2 , the treatment of individuals with clinical signs of the disease u 3 , killing adult vectors with adulticide, u 4 ) are most cost-effective that all the strategies studied in this work.

Conclusion

In this paper, we derived and analysed a model for the control of arboviral diseases with non linear form of infection and complete stage structured model for vectors, and which takes into account a vaccination with waning immunity, treatment, individual protection and vector control strategies (adult vectors, eggs and larvae reduction strategies).

We have begun by calculate the net reproductive number N and the basic reproduction number R 0 , of the basic model (the model without control), and investigate the existence and stability of equilibria. The stability analysis revealed that for N ≤ 1, the trivial equilibrium is globally asymptotically stable. When N > 1 and R 0 < 1, the disease-free equilibrium is locally asymptotically stable. We have found that the model exhibits backward bifurcation. The epidemiological implication of this phenomenon is that for effective eradication and control of diseases, R 0 should be less than a critical values less than one. We have explicitly derived threshold conditions for saddle-node bifurcation in term of the transmission rate, β hv , as well as the basic reproduction number. Then, we have proved, that the disease-induced death is the principal cause of the backward bifurcation phenomenon in model.

Using data from literature related to the transmission dynamics of dengue fever, we also estimated the probability that the model predicts the existence of multiple endemic equilibrium and of the likely stability of the disease-free equilibrium point, through Latin Hypercube Sampling (LHS). The result showed that the model is in an endemic state, since the mean of R 0 is greater than unity. Then, using global sensitivity analysis, we have computed the Partial Rank Correlation Coefficients between R 0 and each parameter of the model. This analysis showed that the basic reproduction number is sensitive to changes in the parameters β vh , the probability of transmission of infection from an infected human to a susceptible vector, β hv , the probability of transmission of infection from an infected vector to a susceptible human, a, the average number of bites, θ, the maturation rate from pupae to adult, µ v , the natural mortality rate of vector, Λ h , the recruitment rate of humans and l, the transfer rate from larvae to pupae, which suggested that the control of the epidemic of arboviral diseases pass through a combination of immunization against arbovirus (vaccination of susceptible humans), individual protection against vector bites, treatment of infected human, vector control through chemical interventions (adulticide and larvicide).

We then considered five time dependent controls as a way out, to ensure the eradication of the disease in a finite time. We performed optimal control analysis of the model. In this light, we addressed the optimal control by deriving and analysing the conditions for optimal eradication of the disease and in a situation where eradication is impossible or of less benefit compared with the cost of intervention, we also derived and analysed the necessary conditions for optimal control of the disease.

From the numerical results and efficiency analysis, as well as, cost-effectiveness analysis, we concluded that the optimal strategy to effectively control arboviral diseases is the combination of vaccination, individual protection, treatment, and other mechanisms of vector control (by chemical intervention-the adulticides). However this conclusion must be taken with caution because of the uncertainties around the parameter values and to the budget/resource limitation. It is also important to note that in most of the work which speak of the optimal control of infectious diseases (see e.g. [START_REF] Adams | Dynamic multidrug therapies for hiv: optimal and sti control approaches[END_REF][START_REF] Buonomo | A simple analysis of vaccination strategies for rubella[END_REF][START_REF] Zaman | Stability analysis and optimal vaccination of an sir epidemic model[END_REF][START_REF] Jung | Optimal control of treatments in a two-strain tuberculosis model[END_REF][START_REF] Yusuf | Optimal control of vaccination and treatment for an sir epidemiological model[END_REF]), and particularly the arboviral diseases [START_REF] Aldila | An optimal control problem arising from a dengue disease transmission model[END_REF][START_REF] Blayneha | Backward bifurcation and optimal control in transmission dynamics of west nile virus[END_REF][START_REF] Moulay | Optimal control of chikungunya disease: larvae reduction,treatment and prevention[END_REF][START_REF] Rodrigues | Vaccination models and optimal control strategies to dengue[END_REF][START_REF] Weverton | A multiobjective optimization approach for combating aedes aegypti using chemical and biological alternated step-size control[END_REF], cost effectiveness analysis, to our knowledge, is not did by the authors. This therefore represents a contribution to the study of optimal control models of arboviral diseases.

In addition, the utilization of a vaccine of small efficacy could have a negative impact on the health of the population. Indeed, for the particular case of dengue, the fact that sequential infections with different strains can cause severe forms of the disease must be taken into account. For instance, it is not currently known if a vaccinated individual, for which the efficacy for a given strain is small, can develop a severe form of the disease when coming into contact with such a strain. Also, its efficiency is higher in children 9-16 years (two thirds are immune) and in individuals who have already been infected. The vaccine appears to contrast against-productive in younger children without the researchers knowing why. The results of clinical trials -which involved more than 40,000 volunteers-were therefore not lifted all the uncertainties about the impact of the vaccine [START_REF] Le Monde Economie | Sanofi lance son vaccin contre la dengue[END_REF]. Therefore, pending the completion of Phase III trials on the efficacy of the vaccine against dengue (Dengvaxia ), and therefore its acceptance by public health organizations such as WHO and the Centre for Disease Control (CDC), it is important to focus on other control mechanisms.

All simulated intervention combinations can be considered cost-effective in the context of available resources for health in countries affected by arboviruses. These results have the potential to help managers control programs against arbovirus infections in high endemicity countries by modifying the implementation of current interventions, or by adding new control mechanisms.

at the Trivial equilibrium is given by

Df (E 0 ) =                 -µ h 0 0 0 0 -aβ hv η v -aβ hv 0 0 0 0 -k 3 0 0 0 aβ hv η v aβ hv 0 0 0 0 γ h -k 4 0 0 0 0 0 0 0 0 0 σ -µ h 0 0 0 0 0 0 0 0 0 0 -µ v 0 0 0 0 θ 0 0 0 0 0 -k 9 0 0 0 0 0 0 0 0 0 γ v -µ v 0 0 0 0 0 0 0 µ b µ b µ b -k 5 0 0 0 0 0 0 0 0 0 s -k 6 0 0 0 0 0 0 0 0 0 l -k 7                 . ( 52
)
The characteristic polynomial of Df (E 0 ) is given by:

P (λ) = (λ -µ h ) 2 (λ -k 3 )(λ -k 4 )(λ -k 9 )(λ -µ v )φ 1 (λ)
where

φ 1 (λ) = λ 4 + A 1 λ 3 + A 2 λ 2 + A 3 λ + A 4 , with A 1 = µ v + k 7 + k 6 + k 5 , A 2 = (k 7 + k 6 + k 5 ) µ v + (k 6 + k 5 ) k 7 + k 5 k 6 , A 3 = ((k 6 + k 5 ) k 7 + k 5 k 6 ) µ v + k 5 k 6 k 7 , A 4 = k 5 k 6 k 7 µ v (1 -N ) .
The roots of P (λ) are

λ 1 = -µ h , λ 1 = -k 3 , λ 2 = -k 4 , λ 3 = -µ v , λ 4 = -k 9
, and the others roots are the roots of φ 1 (λ). Since N < 1, it is clear that all coefficients of φ 1 (λ) are always positive. Now we just have to verify that the Routh-Hurwitz criterion holds for polynomial φ 2 (λ). To this aim, setting

H 1 = A 1 , H 2 = A 1 1 A 3 A 2 , H 3 = A 1 1 0 A 3 A 2 A 1 0 A 4 A 3 , H 4 = A 1 1 0 0 A 3 A 2 A 1 1 0 A 4 A 3 A 2 0 0 0 A 4 = A 4 H 3 .
The Routh-Hurwitz criterion of stability of the trivial equilibrium E 0 is given by

       H 1 > 0 H 2 > 0 H 3 > 0 H 4 > 0 ⇔        H 1 > 0 H 2 > 0 H 3 > 0 A 4 > 0 (53)
We have

H 1 = A 1 > 0, H 2 = A 1 A 2 -A 3 = (k 7 + k 6 + k 5 ) µ 2 v + k 2 7 + (2 k 6 + 2 k 5 ) k 7 + k 2 6 + 2 k 5 k 6 + k 2 5 µ v + (k 6 + k 5 ) k 2 7 + k 6 + k 2 5 2 k 7 + k 5 k 2 6 + k 2 5 k 6 , H 3 = A 1 A 2 A 3 -A 2 1 A 4 -A 2 3 = (k 6 + k 5 ) k 2 7 + k 2 6 + 2k 5 k 6 + k 2 5 k 7 + k 5 k 2 6 + k 2 5 k 6 µ 3 v + µ b lsθ + (k 6 + k ) k 3 7 + 2k 2 6 + 4k 5 k 6 + 2k 2 5 k 2 7 + k 3 6 + 4k 5 k 2 6 + 4k 2 5 k 6 + k 3 5 k 7 + k 5 k 3 6 + 2k 2 5 k 2 6 + k 3 5 k 6 µ 2 v + (2k 7 + 2k 6 + 2k 5 ) µ b lsθ + k 2 6 + 2k 5 k 6 + k 2 5 k 3 7 + k 3 6 + 4k 5 k 2 6 + 4k 2 5 k 6 + k 3 5 k 2 7 + 2k 5 k 3 6 + 4k 2 5 k 2 6 + 2k 3 5 k 6 k 7 + k 2 5 k 3 6 + k 3 5 k 2 6 µv + k 2 7 + (2k 6 + 2k 5 ) k 7 + k 2 6 + 2k 5 k 6 + k 2 5 µ b lsθ + k 5 k 2 6 + k 2 5 k 6 k 3 7 + k 5 k 3 6 + 2k 2 5 k 2 6 + k 3 5 k 6 k 2 7 + k 2 5 k 3 6 + k 3 5 k 2 6 k 7
We always have

H 1 > 0, H 2 > 0, H 3 > 0 and H 4 > 0 if N < 1.
Thus, the trivial equilibrium E 0 is locally asymptotically stable whenever N < 1.

We assume the net reproductive number N > 1. Following the procedure and the notation in [START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF], we may obtain the basic reproduction number R 0 as the dominant eigenvalue of the next-generation matrix [START_REF] Diekmann | Mathematical Epidemiology of Infectious Diseases. Model building, analysis and interpretation[END_REF][START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF]. Observe that model [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF] has four infected populations, namely E h , I h , E v and I v . It follows that the matrices F and V defined in [START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF], which take into account the new infection terms and remaining transfer terms, respectively, are given by

F =       0 0 β hv η v β hv 0 0 0 0 β vh η v N 0 v N 0 h β vh N 0 v N 0 h 0 0 0 0 0 0       and V =     k 3 0 0 0 -γ h k 4 0 0 0 0 k 9 0 0 0 -γ v k 8     .
The dominant eigenvalue of the next-generation matrix F V -1 is given by [START_REF] Cannon | An epidemiology model suggested by yellow fever[END_REF]. The local stability of the disease-free equilibrium E 1 is a direct consequence of Theorem 2 in [START_REF] Van Den Driessche | Reproduction numbers and the sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF]. This ends the proof.

B Proof of Theorem 2

Setting Y = X-T E with X = (S h , E h , I h , R h , S v , E v , I v , E, L, P ) T , A 88 = k 5 + µ b S v + E v + I v K E ,
and

A 99 = k 6 + s E Γ L
. we can rewrite (3) in the following manner

dY dt = B(Y )Y (54)
where

B(Y ) =                    -(λ h + µ h ) 0 0 0 0 - aβ hv η v S 0 h N h - aβ hv S 0 h N h 0 0 0 λ h -k 3 0 0 0 aβ hv η v S 0 h N h aβ hv S 0 h N h 0 0 0 0 γ h -k 4 0 0 0 0 0 0 0 0 0 σ -µ h 0 0 0 0 0 0 0 0 0 0 -(λ v + µ v ) 0 0 0 0 θ 0 0 0 0 λ v -k 9 0 0 0 0 0 0 0 0 0 γ v -µ v 0 0 0 0 0 0 0 µ b µ b µ b -A 88 0 0 0 0 0 0 0 0 0 s -A 99 0 0 0 0 0 0 0 0 0 l -k 7                    .
It is clear that Y = (0, 0, 0, 0, 0, 0, 0, 0, 0, 0) is the only equilibrium. Then it suffices to consider the following Lyapunov function L(Y ) =< g, Y > were g = 1, 1, 1, 1, 1, 1, 1,

k 8 µ b , k 5 k 8 µ b s , k 5 k 6 k 8 µ b sl . Straightforward computations lead that L(Y ) =< g, Ẏ > def = < g, B(Y )Y > = -µ h Y 1 -µ h Y 2 -(µ h + δ)Y 3 -µ h Y 4 - k 8 K E (Y 5 + Y 6 + Y 7 ) - k 5 k 8 µ b K L Y 8 Y 9 + θ 1 - 1 N Y 10 .
We have L(Y ) < 0 if N ≤ 1 and L(Y ) = 0 if Y i = 0, i = 1, 2, . . . , 10 (i.e S h = S 0 h and

E h = I h = R h = S v = E v = I v = E = L = P = 0)
. Moreover, the maximal invariant set contained in L| L(Y ) = 0 is {(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)}. Thus, from Lyapunov theory, we deduce that {(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)} and thus, E 0 , is GAS if and only if N ≤ 1.

C Proof of Theorem 3

Proof. In order to determine the existence of endemic equilibria, i.e., equilibria with all positive components, say

E * * = (S * h , V * h , E * h , I * h , R * h , S * v , E * v , I *
v , E, L, P ) , we have to look for the solution of the algebraic system of equations obtained by equating the right sides of system (6) to zero.

Solving the equations in the model ( 6) in terms of λ * h and λ * v , gives

S * h = Λ h µ h + λ * h , E * h = λ * h S * h k 3 , I * h = γ h λ * h S * h k 3 k 4 , R * h = σγ h λ * h S * h µ h k 3 k 4 , (55)
and

S * v = θP (λ * v + k 8 ) , E * v = θP λ * v k 9 (λ * v + k 8 ) , I * v = γ v θP λ * v k 8 k 9 (λ * v + k 8 ) , E = µ b θΓ E P (k 5 k 8 Γ E + µ b θP ) , L = µ b θsΓ E Γ L P k 6 Γ L (k 5 k 8 Γ E + µ b θP ) + sµ b θΓ E P , ( 56
)
where P is solution of the following equation

f (P ) = -k 7 P [µ b θ(sΓ E + k 6 Γ L )P -k 5 k 6 k 8 Γ E Γ L (N -1)] = 0 (57)
A direct resolution of the above equation give P = 0 or

P = k 5 k 6 k 8 Γ E Γ L (N -1) µ b θ(sΓ E + k 6 Γ L ) .
Note that P = 0 corresponds to the trivial equilibrium E 0 . Now we consider P > 0 i.e. N > 1. Substituting ( 55) and ( 56) into the expression of λ * h and λ * v and simplifying, lead the non-zero equilibria of the basic model ( 6) satisfy the quadratic equation

k 9 µ b Λ h (sΓ E + k 6 Γ L ) d 2 (λ * h ) 2 + d 1 λ * h + d 0 = 0, (58)
where d i , i = 0, 1, 2, are given by

d 2 = -k 2 [k 10 aµ h β vh + k 2 k 8 ] < 0, (59a) d 1 = k 2 3 k 2 4 k 8 µ h (R 2 0 -R 2 c ), (59b)
d 0 = k 2 3 k 2 4 k 8 µ 2 h R 2 0 -1 . ( 59c
)
and R 0 and R c are given by ( 10) and [START_REF] Carr | Applications of Centre Manifold Theory[END_REF], respectively. This equation may be simply analyzed through the Descartes' rule of signs. First of all, note that d 2 is negative. Therefore the following cases are possible:

1. There is a unique endemic equilibrium if d 0 > 0;

4. There are no endemic equilibria otherwise.

We observe d 0 > 0 is equivalent to R 0 > 1 so statement (i) of Theorem 3 is equivalent to point (a).

When R 0 = 1, d 0 = 0. We observe that d 1 > 0 is equivalent to R c < R 0 . Therefore, when R 0 = 1 and R c < 1, d 0 = 0 and d 1 > 0, so statement (ii) a) of Theorem 3 follows from statement (b) above. Since the condition d 0 = 0 does not appear elsewhere in statements (a), (b), or (c) above, statement (ii) b) of Theorem 3 follows from statement (d ) above.

When R 0 < 1, d 0 < 0, and when R c < R 0 , d 1 > 0. We also note that for R 0 < 1, when

d 1 > 0, ψ ≤ 0 because ψ > 0 is equivalent to d 1 < 0. Indeed, ψ > 0 ⇐⇒ k 10 aµ h β vh > δγ h k 8 ⇐⇒ k 10 aµ h β vh k 9 + 2k 8 k 9 k 2 > δγ h k 8 k 9 + 2k 8 k 9 k 2 ⇐⇒ a 2 β hv β vh k 3 k 4 k 10 k 11 µ h N 0 v -k 3 k 4 µ h N 0 h (k 10 aµ h β vh k 9 + 2k 8 k 9 k 2 ) < k 3 k 4 k 8 k 9 µ h N 0 h (R 2 0 -1)k 3 k 4 -k 3 k 4 + δγ h ) ⇐⇒ a 2 β hv β vh k 3 k 4 k 10 k 11 µ h N 0 v -k 3 k 4 µ h N 0 h (k 10 aµ h β vh k 9 + 2k 8 k 9 k 2 ) < k 3 k 4 k 8 k 9 µ h N 0 h (R 2 0 -1)k 3 k 4 -k 2 ) ⇐⇒ d 1 < k 3 k 4 k 8 µ h (R 2 0 -1)k 3 k 4 -k 2 ) < 0, since R 0 < 1. (62)
Consequently, we show that d 2 1 -4d 2 d 0 = 0 is equivalent to,

ρ 2 R 4 0 + ρ 1 R 2 0 + ρ 0 = 0, (63)
where

ρ 2 = k 4 3 k 4 4 k 2 8 µ 2 h , (64a) ρ 1 = 2k 2 3 k 2 4 k 8 µ 2 h [k 2 (k 10 aµ h β vh -k 8 δγ h ) -(k 10 aµ h β vh + k 8 k 2 )δγ h ] , (64b)
ρ 0 = k 2 3 k 2 4 k 2 10 a 2 µ 4 h β 2 vh . (64c)
We again use Descartes' rule of signs to analyse equation [START_REF] Parks | Planning social mobilization and communication for dengue fever prevention and control[END_REF]. The discriminant of ( 63) is ∆ r = ρ 2 1 -4ρ 2 ρ 0 , and can be written

∆ r = -16k 4 3 k 4 4 k 2 8 k 2 µ 4 h δγ h (k 10 aµ h β vh + k 8 k 2 ) (k 10 aµ h β vh -k 8 δγ h )
Since ρ 2 > 0 and ρ 0 > 0, equation ( 63) allows real positive solutions if and only if ρ 1 < 0 and ∆ r ≥ 0. Now, we express the obtained inequalities in terms of the quantities ( 11)- [START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF]. To this aim, we note that ∆ r ≥ 0 is equivalent to ψ ≤ 0. From the definition of ψ (11) and ρ 1 (64b), this implies that ρ 1 < 0. Therefore, equation ( 63) has exactly two positive solutions given by ( 13) and [START_REF] Castillo-Chavez | Dynamical models of tuberculosis and their applications[END_REF]. Therefore, statement (iii) b) follows from statement (b) above.

Similarly, d 2 1 -4d 2 d 0 > 0, with d 1 > 0 and d 0 < 0 written in terms of the basic reproduction number, is equivalent to

R 0 < R 1b or R 0 > R 2b , (65)
so statement (iii) a) follows from statement (c) above. Finally, statement (iii) c) then follows from statement (d ) above. Thus Theorem 3 is established.

D Derivation of formula (21)

The local bifurcation analysis near the bifurcation point (β hv = β * hv ) is then determined by the signs of two associated constants, denoted by A 1 and A 2 , defined by

A 1 = 10 k,i,j=1 v k w i w j ∂ 2 f k (0, 0) ∂x i ∂x j and A 2 = 10 k,i=1 v k w i ∂ 2 f k (0, 0) ∂x i ∂φ (66)
with φ = β hvβ * hv . It is important to note that in f k (0, 0), the first zero corresponds to the disease-free equilibrium, E 1 , for the system [START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF]. Since β hv = β * hv is the bifurcation parameter, it follows from φ = β hvβ * hv that φ = 0 when β hv = β * hv which is the second component in f k (0, 0). Using Eqs. ( 17) and [START_REF] Cruz-Pacheco | Seasonality and outbreaks in west nile virus infection[END_REF] in Eq. ( 19), we obtain

A 1 = v 2 10 i,j=1 w i w j ∂ 2 f 2 (0, 0) ∂x i ∂x j + v 6 10 i,j=1 w i w j ∂ 2 f 6 (0, 0) ∂x i ∂x j and A 2 = v 2 10 i=1 w i ∂ 2 f 2 (0, 0) ∂x i ∂φ . (67)
Let A

= 10 i,j=1 w i w j ∂ 2 f 2 (0, 0) ∂x i ∂x j = w 1 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 1 ∂x j + w 2 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 2 ∂x j + w 3 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 3 ∂x j + w 4 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 4 ∂x j + w 5 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 5 ∂x j + w 6 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 6 ∂x j + w 7 10 i,j=1 w j ∂ 2 f 2 (0, 0) ∂x 7 ∂x j = w 2 w 6 ∂ 2 f 2 ∂E h ∂E v (0, 0) + w 7 ∂ 2 f 2 ∂E h ∂I v (0, 0) + w 3 w 6 ∂ 2 f 2 ∂I h ∂E v (0, 0) + w 7 ∂ 2 f 2 ∂I h ∂I v (0, 0) + w 4 w 6 ∂ 2 f 2 ∂R h ∂E v (0, 0) + w 7 ∂ 2 f 2 ∂R h ∂I v (0, 0) + w 6 w 2 ∂ 2 f 2 ∂E v ∂E h + w 3 ∂ 2 f 2 ∂E v ∂I h + w 4 ∂ 2 f 2 ∂E v ∂R h + w 7 w 2 ∂ 2 f 2 ∂I v ∂E h + w 3 ∂ 2 f 2 ∂I v ∂I h + w 4 ∂ 2 f 2 ∂I v ∂R h = w 2 - aβ * hv η v N 0 h w 6 - aβ * hv N 0 h w 7 + w 3 - aβ * hv η v N 0 h w 6 - aβ * hv N 0 h w 7 + w 4 - aβ * hv η v N 0 h w 6 - aβ * hv N 0 h w 7 + w 6 - aβ * hv η v N 0 h w 2 - aβ * hv η v N 0 h w 3 - aβ * hv η v N 0 h w 4 + w 7 - aβ * hv N 0 h w 2 - aβ * hv N 0 h w 3 - aβ * hv N 0 h w 4 = - aβ * hv N 0 h (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ) - aβ * hv N 0 h (w 2 + w 3 + w 4 ) (η v w 6 + w 7 ) = -2 aβ * hv N 0 h (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ), (1) 1
A (2) 1 = 10 i,j=1 w i w j ∂ 2 f 6 (0, 0) ∂x 1 ∂x j = w 1 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 1 ∂x j + w 2 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 2 ∂x j + w 3 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 3 ∂x j + w 4 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 4 ∂x j + w 5 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 5 ∂x j + w 6 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 6 ∂x j + w 7 10 i,j=1 w j ∂ 2 f 6 (0, 0) ∂x 7 ∂x j = w 1 w 2 ∂f 6 ∂S h ∂E h (0, 0) + w 3 ∂f 6 ∂S h ∂I h (0, 0) + w 2 w 1 ∂f 6 ∂E h ∂S h (0, 0) + w 2 ∂f 6 ∂E 2 h (0, 0) + w 3 ∂f 6 ∂E h ∂I h (0, 0) + w 4 ∂f 6 ∂E h ∂R h (0, 0) + w 5 ∂f 6 ∂E h ∂Sv (0, 0) + w 3 w 1 ∂f 6 ∂I h ∂S h (0, 0) + w 2 ∂f 6 ∂I h ∂E h (0, 0) + w 3 ∂f 6 ∂I 2 h (0, 0) + w 4 ∂f 6 ∂I h ∂R h (0, 0) + w 5 ∂f 6 ∂I h ∂Sv (0, 0) + w 4 w 2 ∂f 6 ∂R h ∂E h (0, 0) + w 3 ∂f 6 ∂R h ∂I h (0, 0) + w 5 w 2 ∂f 6 ∂Sv∂E h (0, 0) + w 3 ∂f 6 ∂Sv∂I h (0, 0) = w 1 - aβ vh η h S 0 v (N 0 h ) 2 w 2 - aβ vh S 0 v (N 0 h ) 2 w 3 + w 2 - aβ vh η h S 0 v (N 0 h ) 2 w 1 - and
aβ vh η h S 0 v (N 0 h ) 2 w 2 - aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 3 - aβ vh η h S 0 v (N 0 h ) 2 w 4 + aβ vh η h N 0 h w 5 + w 3 - aβ vh S 0 v (N 0 h ) 2 w 1 - aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 -2 aβ vh η h S 0 v (N 0 h ) 2 w 3 - aβ vh S 0 v (N 0 h ) 2 w 4 + aβ vh N 0 h w 5 + w 4 - aβ vh η h S 0 v (N 0 h ) 2 w 2 - aβ vh S 0 v (N 0 h ) 2 w 3 + w 5 aβ vh η h N 0 h w 2 + aβ vh N 0 h w 3 = - aβ vh η h S 0 v (N 0 h ) 2 w 1 w 2 - aβ vh S 0 v (N 0 h ) 2 w 1 w 3 - aβ vh η h S 0 v (N 0 h ) 2 w 1 w 2 -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 2 - aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 w 3 - aβ vh η h S 0 v (N 0 h ) 2 w 2 w 4 + aβ vh η h N 0 h w 2 w 5 - aβ vh S 0 v (N 0 h ) 2 w 1 w 3 - aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 w 3 - 2
aβ vh η h S 0 v (N 0 h ) 2 w 2 3 - aβ vh S 0 v (N 0 h ) 2 w 3 w 4 + aβ vh N 0 h w 3 w 5 - aβ vh η h S 0 v (N 0 h ) 2 w 2 w 4 - aβ vh S 0 v (N 0 h ) 2 w 3 w 4 + aβ vh η h N 0 h w 2 w 5 + aβ vh N 0 h w 3 w 5 = -2 aβ vh η h S 0 v (N 0 h ) 2 w 1 w 2 -2 aβ vh S 0 v (N 0 h ) 2 w 1 w 3 -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 2 -2 aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 w 3 -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 w 4 -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 3 - 2
aβ vh S 0 v (N 0 h ) 2 w 3 w 4 + 2 aβ vh η h N 0 h w 2 w 5 + 2 aβ vh N 0 h w 3 w 5 = -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 2 -2 aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 w 3 -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 w 4 - 2
aβ vh η h S 0 v (N 0 h ) 2 w 2 3 - 2
aβ vh S 0 v (N 0 h ) 2 w 3 w 4 + 2 aβ vh N 0 h η h w 2 w 5 + 2 aβ vh N 0 h w 3 w 5 -2 aβ vh S 0 v (N 0 h ) 2 η h w 1 w 2 -2 aβ vh S 0 v (N 0 h ) 2 w 1 w 3 = -2 aβ vh η h S 0 v (N 0 h ) 2 w 2 2 - 2
aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 w 3 - 2
aβ vh η h S 0 v (N 0 h ) 2 w 2 w 4 - 2
aβ vh η h S 0 v (N 0 h ) 2 w 2 3 - 2
aβ vh S 0 v (N 0 h ) 2 w 3 w 4 + 2 aβ vh N 0 h (η h w 2 + w 3 ) w 5 -2 aβ vh S 0 v (N 0 h ) 2 η h w 1 w 2 -2 aβ vh S 0 v (N 0 h ) 2 w 1 w 3 = -2 aβ vh η h S 0 v (N 0 h ) 2 w 2
According to [START_REF] Cruz-Pacheco | Modelling the dynamics of west nile virus[END_REF] and ( 18), we have ζ 1 > 0 and ζ 2 > 0.

We then have

A 2 = v 2 10 i=1 w i ∂ 2 f 2 (0, 0) ∂x i ∂φ = v 2 w 6 ∂ 2 f 2 ∂E v ∂φ (0, 0) + w 7 ∂ 2 f 2 ∂E v ∂φ (0, 0) = v 2 η h w 6 aS 0 v N 0 h + w 7 aS 0 v N 0 h = aS 0 v N 0 h (η h w 6 + w 7 ) v 2 .

E Proof of Theorem 5

We follow the approach given in [START_REF] Safan | Vaccination based control of infections in sirs models with reinfection: special reference to pertussis[END_REF]. At this aim, note that equation ( 58) may be written as

F (β hv , λ h ) := d 2 (λ * h ) 2 + d 1 λ * h + d 0 = 0, (70)
where d 2 , d 1 and d 0 are the same coefficients as in [START_REF] Stein | Large sample properties of simulations using latin hypercube sampling[END_REF]. Thus, the positive endemic equilibria of model ( 6) are obtained by solving (70) for positive λ * h and substituting the results into [START_REF] Scott | Vector dynamics and transmission of dengue virus: implications for dengue surveillance and prevention strategies: vector dynamics and dengue prevention[END_REF]. Clearly, the coefficient d 2 , of (70), is always negative while d 1 and d 0 may change sign. Therefore, there is a single endemic equilibrium if and only if d 0 > 0, which correspond to R 0 > 1. There are two endemic equilibria if and only if d 0 < 0, d 1 > 0 and d 2 1 -4d 2 d 0 > 0. Now, first remember that d 0 < 0 (i. e. R 0 < 1) is equivalent to β hv < β * hv , where β * hv is given at Eq. ( 16).

Then, inequality d 1 > 0, is equivalent to

β hv > β ( 71
)
where β is given by [START_REF] Dumont | Vector control for the chikungunya disease[END_REF]. Finally, equation d 2 1 -4d 2 d 0 = 0, in terms of β hv , is equivalent to

α 2 β 2 hv + α 1 β hv + α 0 = 0, ( 72
)
where α 2 = k Under condition (73), we conclude that α 1 < 0. Thus, Eq. (72) admits exactly two positive solutions which are given by

β ± = -α 1 ± √ ∆ 2α 2 = N 0 h k 3 k 4 k 10 k 11 a 2 N 0 v β vh δγ h (aµ h β vh k 10 + k 2 k 8 ) ± (-k 2 ψ) 2
Thus, condition d 2 1 -4d 2 d 0 > 0 written in the terms of β hv is equivalent to

β hv < β -or β hv > β + .
and the inequalities (25) then follow.

Remark 3. Note that condition (25) is equivalent to condition [START_REF]Dengue and severe dengue[END_REF]. Therefore, Theorem 5 is equivalent to Theorem 3.

F Proof of Theorem 6

Considering the model ( 6) without disease-induced death in human, and applying the same procedure as appendix 3, we obtain that the non-zero equilibria of the basic model ( 6) satisfy the linear equation (sΓ E + k 6 Γ L )(p 1 λ * h + p 0 ) = 0, where p 1 = µ b Λ h k 9 (k 2 aµ h β vh + k 3 k 8 (µ h + σ)) and p 0 = -µ h k 3 k 4 k 8 k 9 µ b Λ h R 2 0,δ=0 -1 . Clearly, p 1 > 0 and p 0 ≥ 0 whenever R 0,δ=0 ≤ 1, so that λ * h = -p 0 p 1 ≤ 0. Therefore, the model ( 6) without disease-induced death in human, has no endemic equilibrium whenever R 0,δ=0 ≤ 1. The above result suggests the impossibility of backward bifurcation in the model ( 6) without disease-induced death, since no endemic equilibrium exists when R 0,δ=0 < 1 (and backward bifurcation requires the presence of at least two endemic equilibria when R 0,δ=0 < 1) [START_REF] Garba | Backward bifurcations in dengue transmission dynamics[END_REF][START_REF] Sharomi | Role of incidence function in vaccine-induced backward bifurcation in some hiv models[END_REF].

To completely rule out backward bifurcation in model ( 6), we use the direct Lyapunov method to prove the global stability of the DFE. Define the positively-invariant and attracting region

D 1 = (S h , E h , I h , R h , S v , E v , I v , E, L, P ) ∈ D : S h ≤ N 0 h ; S v ≤ N 0 v .
Consider the Lyapunov function

G = q 1 E h + q 2 I h + q 3 E v + q 4 I v ,
where

q 1 = 1 k 3 ; q 3 = τ 1 S 0 h k 3 k 8 k 11 k 9 , q 2 = τ 1 S 0 h k 3 k 8 k 11 ζ 2 S 0 v k 4 k 9 , q 4 = τ 1 S 0 h k 3 k 8 .
and we have set τ 1 = µ h β hv Λ h and τ 2 = µ h β vh Λ h . The derivative of G is given by Ġ = q 1 Ėh + q 2 İ h + q 3 Ėv + q 4 İv = q 1 (λ h S hk 3 E h ) + q 2 (γ h E hk 4 I h ) + q 3 (λ v S vk 9 E v ) + q 4 (γ v E vk 8 I v ) = q 1 τ 1 S h (η v E v + I v )q 3 k 9 E v + q 4 γ v E vq 4 k 8 I v + q 3 τ 2 S v (η h E h + I h )q 1 k 3 E h + q 2 γ h E hq 2 k 4 I h = (q 1 τ 1 S h η v + q 4 γ vq 3 k 9 )E v + (q 1 τ 1 S hq 4 k 8 )I v + (q 3 τ 2 S v η h + q 2 γ hq 1 k 3 )E h + (q 3 τ 2 S vq 2 k 4 )I h ≤ (q 1 τ 1 S 0 h η v + q 4 γ vq 3 k 9 )E v + (q 1 τ 1 S 0 hq 4 k 8 )I v + (q 3 τ 2 S 0 v η h + q 2 γ hq 1 k 3 )E h + (q 3 ζ 2 S 0 vq 2 k 4 )I h , since S h ≤ S 0 h , S v ≤ S 0 v in D 1 .

Replacing q i , i = 1, . . . , 4, by their value gives after straightforward simplifications Ġ ≤ R 2 0,δ=0 -1 E h

We have Ġ ≤ 0 if R 0,δ=0 ≤ 1, with Ġ = 0 if R 0,δ=0 = 1 or E h = 0. Whenever E h = 0, we also have I h = 0, E v = 0 and I v = 0. Substituting E h = I h = E v = I v = 0 in the first, fourth and fifth equation of Eq. ( 6) with δ = 0 gives S h (t) → S 0 h = N 0 h , R h (t) → 0, and S v (t) → S 0 v = N 0 v as t → ∞. Thus

[S h (t), E h (t), I h (t), R h (t), S v (t), E v (t), I v (t), E(t), L(t), P (t)] → (N 0 h , 0, 0, 0, N 0 v , 0, 0, E, L, P ) as t → ∞.

It follows from the LaSalle's invariance principle [START_REF] Hale | Ordinary Differential Equations[END_REF][START_REF] Lasalle | Stability theory for ordinary differential equations[END_REF][START_REF]The stability of dynamical systems[END_REF] that every solution of (6) (when R 0,δ=0 ≤ 1), with initial conditions in D 1 converges to E 1 , as t → ∞. Hence, the DFE, E 1 , of model ( 6) without disease-induced death, is GAS in D 1 if R 0,δ=0 ≤ 1.

Figure 1 :

 1 Figure 1: Schematic of the vector-borne epidemic model with development stage of vectors.

Theorem 4 .

 4 The basic model (6) exhibits a backward bifurcation at R 0 = 1 whenever A 1 > 0 (i.e., ζ 1 > ζ 2 ). If the reversed inequality holds, then the bifurcation at R 0 = 1 is forward.

Figure 2 :

 2 Figure 2: The backward bifurcation curves for model system (6) in the (R 0 , I * h ) and (R 0 , I * v ) planes. The parameter β hv varied in the range [0, 0.0877] to allow R 0 to vary in the range [0, 1.5]. Two endemic equilibrium points coexist for values of R 0 in the range (0.2671, 1)(corresponding to the range (0.0028, 0.0390) of β hv ). The notation EE and DFE stand for endemic equilibrium and disease-free equilibrium, respectively. Solid lines represent stable equilibria and dash lines stand for unstable equilibria.

Figure 3 :

 3 Figure 3: Solutions of model (6) of the number of infectious humans, I h , and the number of infectious vectors, I v , for parameter values given in the bifurcation diagram in Figure2with β hv = 0.008, so R 0 = 0.4359 < 1, for two different set of initial conditions. The first set of initial conditions (corresponding to the dotted trajectory) is S h = 700, E h = 220, I h = 15, R h = 60, S v = 3000, E v = 400, I v = 120, E = 10000, L = 5000 and P = 3000. The second set of initial conditions (corresponding to the solid trajectory) is S h = 733650, E h = 220, I h = 15, R h = 60, S v = 3000, E v = 400, I v = 120, E = 10000, L = 5000 and P = 3000. The solution for initial condition 1 approaches the locally asymptotically stable DFE point, while the solution for initial condition 2 approaches the locally asymptotically stable endemic equilibrium .

FreqenciesFigure 4 :

 4 Figure 4: Sampling distribution of R 0 from 5,000 runs of Latin hypercube sampling. The mean of R 0 is 1.9583 and the standard deviation is 1.8439. Furthermore, P(R 0 ≥ 1) = 64.48%.

Figure 5 :

 5 Figure 5: Partial rank correlation coefficients for R 0

Figure 6 :

 6 Figure 6: Control functions.

Figure 7 :

 7 Figure 7: Simulation results of optimal control model[START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with individual protection, treatment and adulticide (u 1 = 0, u 2 = 0, u 3 = 0, u 4 = 0 ).

Figure 8 :

 8 Figure 8: Simulation results of optimal control model[START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with individual protection, treatment and larvicide (u 1 = 0, u 2 = 0, u 3 = 0, u 5 = 0 ).

Figure 9 :

 9 Figure9: Simulation results of optimal control model[START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with treatment, adulticide and larvicide (u 1 = 0, u 3 = 0, u 4 = 0, u 5 = 0 ).

Figure 10 :

 10 Figure 10: Simulation results of optimal control model[START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using optimal vaccination combined with individual protection, adulticide and larvicide (u 1 = 0, u 2 = 0, u 4 = 0, u 5 = 0 ).

Figure 11 :

 11 Figure 11: Simulation results of optimal control model[START_REF] Fleming | Deterministic and Stochastic Optimal Control[END_REF] showing the effect of using the combination of all the five controls (u i = 0, i = 1, . . . , 5).

2 . 60 ) 3 .

 2603 There is a unique endemic equilibrium if (d 1 > 0 and d 0 = 0) or (d 1 > 0 and d 0 < 0 and d 2 1 -4d 2 d 0 = 0); (There are two endemic equilibria if d 1 > 0 and d 0 < 0 and d 2 1 -4d 2 d 0 > 0;

Table 1 :

 1 The state variables of model[START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF].

	Humans	Aquatic Vectors	Adult Vectors
	S h : Susceptible	E: Eggs	S

v : Susceptible E h : Infected in latent stage L: Larvae E v Infected in latent stage I h : Infectious P : Pupae I v Infectious R h : Resistant (immune)

Table 2 :

 2 Description and baseline values/range of parameters of model[START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF]. The baseline values refer to dengue fever transmission.

	Parameter Description	Baseline	Sources
			value/range
	Λ h	Recruitment rate of humans	2.5 day -1	[29]
	µ h	Natural mortality rate in humans	1 (67×365) day -1	[29]
	a	Average number of bites	1 day -1	[4, 29]
	β hv	Probability of transmission of	0.1, 0.75 day -1 [4, 29]
		infection from an infected vector
		to a susceptible human
	γ h δ	Progression rate from E h to I h Disease-induced death rate	1 15 , 1 3 10 -3 day -1 day -1	[23, 55] [29]
	σ	Recovery rate for humans	0.1428 day -1	[4, 29]
	η h ,η v	Modifications parameter	[0, 1)	[29]
	µ v γ v β vh	Natural mortality rate of vectors Progression rate from E v to I v Probability of transmission of	1 30 , 1 14 1 21 , 1 2 day -1 day -1 0.1, 0.75 day -1 [4, 29] [4, 29] [23, 55]
		infection from an infected human
		to a susceptible vector
	θ	Maturation rate from pupae to adult 0.08 day -1	[23, 43, 44]
	µ b	Number of eggs at each deposit	6 day -1	[23, 43, 44]
	Γ E	Carrying capacity for eggs	10 3 , 10 6	[4, 43]
	Γ L µ E	Carrying capacity for larvae Eggs death rate	5 × 10 2 , 5 × 10 5 [4, 43] 0.2 or 0.4 [44]
	µ L	Larvae death rate	0.2 or 0.4	[44]
	µ P	Pupae death rate	0.4	Assumed
	s	Transfer rate from eggs to larvae	0.7 day -1	[44]
	l	Transfer rate from larvae to pupae	0.5 day -1	[43]
	Remark 1. (i) It is important to note that, in the case of other arboviral diseases (e.g, Chikun-
	gunya), the exposed humans and vectors do not play any role in the infectious process, in this
	case η h

Table 3 :

 3 Partial Rank Correlation Coefficients between R 0 and each parameters of model[START_REF] Antonio | Optimal and sub-optimal control in dengue epidemics[END_REF].

	Parameter Correlation Parameter Correlation Parameter Correlation
		Coefficients		Coefficients		Coefficients
	β vh	0.7345	Γ L	0.3813	µ h	0.1717
	β hv	0.7285	Γ E	0.3698	γ v	0.0872
	a	0.6454	s	0.3733	η v	0.06891
	θ	0.6187	µ P	-0.2565	µ L	-0.0556
	µ v	-0.5521	η h	0.2331	µ b	0.0531
	Λ h	-0.5435	γ h	-0.2204	µ E	-0.0022
	l	0.5144	σ	-0.1867	δ	0.0003

Table 4 :

 4 Value of parameters using in numerical simulations.

	Parameter value Parameter value Parameter value Parameter value
	µ v a Λ h	1 30 1 2.5	l µ E µ b	0.5 0.2 6	α 2 µ h θ	0.5 1 67 * 365 0.08	γ h γ v µ P	1 14 1 21 0.4
	β hv	0.75			σ	0.1428	η v	0.35
	β vh	0.75	ω	0.05	µ L	0.4	δ	10 -3
	Γ E	10000	s	0.7	η 1	0.001	η 2	0.3
	Γ L	5000	η h	0.35	c m	0.2	α 1	0.5

Table 5 :

 5 Numerical values for the cost functional parameters.

	Parameters Value Source	Parameters Value	Source
	D 1 :	10,000 [44]	B 1	10 Assumed
	D 2 :	10,000 [44]	B 2 :	10	[44]
	D 3 :	5000	Assumed	B 3 :	10	[44]
	D 4 :	1	[44]	B 4 :	10 Assumed
				B 5	10	[44]

Table 6 :

 6 Initial conditions.

	Human states Initial value Adult Vector Initial value Aquatic states Initial value
			states
	S h 0 :	700	S v 0	3000	E 0	10000
	E h 0 :	220	E v 0	400	L 0	5000
	I h 0 :	100	I v 0	120	P 0	3000
	R h 0 :	60

Table 7 :

 7 Description of the different Control strategies.

	Strategy Description
	Z 1	Vaccine combined with individual protection, treatment and adulticide
	Z 2

table 8

 8

			3 0%	Z 1	490.6350	88.048451%
	Z 2	528.9018	87.116296% Z 3	580.5772	85.857517%
	Z 4	1.8391 × 10 3 55.200721% Z	490.6350	88.048451%
	From

, we can conclude that the best strategies are Z 1 and Z.

Remark 2. The same reasoning (about efficiency analysis) can be done for vector population, by replacing I h by N v , and A h by A Nv in Eq. (

46

) and (47), respectively.

Table 9 :

 9 Cost Effectiveness of different strategies.

	Strategies Total number of infected Total infection Total cost($)
		individuals	averted
	Z 4 Z 3 Z 2 Z 1 Z	1839 581 529 491 491	87538 88796 88848 88886 88886	7.6218 × 10 8 7.6093 × 10 8 1.6452 × 10 9 7.6081 × 10 8 7.6081 × 10 8

Table 10 :

 10 Strategy Z 4 vs Strategy Z 3 .

	Strategies Total infection averted Total cost ($)
	Z 4 Z 3	87538 88796	7.6218 × 10 8 7.6093 × 10 8
	ICER(Z 4 ) = ICER(Z 3 ) =	7.6218 × 10 8 87538 (7.6093 -7.6218) × 10 8 = 8706.8473 88796 -87538 = -993.6407

Table 11 :

 11 Strategy Z 3 vs Strategy Z 2 .

	Strategies Total infection averted Total cost ($)
	Z 3 Z 2	88796 88848	7.6093 × 10 8 1.6452 × 10 9

Table 12 :

 12 Strategy Z 3 vs Strategy Z 1 .

	Strategies Total infection averted Total cost ($)
	Z 3 Z 1	88796 88886	7.6093 × 10 8 7.6081 × 10 8
	This leads to the following values for the ICER,
	ICER(Z 3 ) = ICER(Z 1 ) =	7.6093 × 10 8 88796 7.6081 × 10 8 -7.6093 × 10 8 = 8569.4175 88886 -88796	= -1333.3333

Table 13 :

 13 Strategy Z 1 vs Strategy Z.

	Strategies Total infection averted Total cost ($)
	Z 1 Z	88886 88886	7.6081 × 10 8 7.6081 × 10 8

  It follows then, after some algebraic computations, thatA 1 = v 2 -2 aβ * (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ) + v 6 -2 aβ vh S 0 v (N 0 h ) 2 η h w 2 2 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 -2 aβ vh (k 8 + γ v ) γ v N 0 h (η h w 2 + w 3 ) w 7 +2 k 7 K 2 K 4 lK 1 K 3 (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ) + v 6 -2 aβ vh S 0 v (N 0 h ) 2 η h w 2 2 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 -2aβ vh (k 8 + γ v ) γ v N 0 h (η h w 2 + w 3 ) w 7 + v 6 2 k 7 K 2 K 4 lK 1 K 3 (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 = v 6 2 k 7 K 2 K 4 lK 1 K 3 (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ) + v 6 -2 aβ vh S 0 v (N 0 h ) 2 η h w 2 2 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 -2aβ vh (k 8 + γ v ) γ v N 0 h (η h w 2 + w 3 ) w 7 = v 6 2 k 7 K 2 K 4 lK 1 K 3 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 + 2 aβ vh (k 8 + γ v ) γ v N 0 (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 , ζ 2 = v 2 2 aβ * (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ) + v 6 2 aβ vh S 0 v (N0h ) 2 η h w 2 2 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 2 aβ vh (k 8 + γ v ) γ v N 0 h (η h w 2 + w 3 ) w 7 .

				hv
			N 0 h
							aβ vh
								N 0 h
	= v 2 -2	aβ * hv N 0 h
						aβ vh
						N 0 h
						aβ vh
						N 0 h
	+ v 2 -2	aβ * hv N 0 h
						aβ vh
						N 0
	-v 6 2	aβ vh S 0 v (N 0 h ) 2 η h w 2 2 (68)
	where we have set
		ζ 1 = v 6 2	k 7 K 2 K 4 lK 1 K 3	aβ vh N 0 h
								hv
								N 0 h
								2 2 -2	aβ vh S 0 v (N 0 h ) 2 (η h + 1)w 2 w 3 -2	aβ vh η h S 0 v (N 0 h ) 2 w 2 w 4 -2	aβ vh η h S 0 v (N 0 h ) 2 w 2 3 -2	aβ vh S 0 v (N 0 h ) 2 w 3 w 4
		+ 2	aβ vh N 0 h	(η h w 2 + w 3 ) -	k 8 + γv γv	w 7 +	k 7 K 2 K 4 lK 1 K 3	w 10 -2	aβ vh S 0 v (N 0 h ) 2 η h w 1 w 2 -2	aβ vh S 0 v (N 0 h ) 2 w 1 w 3
		= -2	aβ vh S 0
		+ 2	k 7 K 2 K 4 lK 1 K 3	aβ vh N 0

v (N 0 h ) 2 η h w 2 2 + (η h + 1)w 2 w 3 + η h w 2 w 4 + η h w 2 3 + w 3 w 4 -

2

aβ vh (k 8 + γv) γvN 0 h (η h w 2 + w 3 ) w 7 h (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 h (η h w 2 + w 3 ) w 10 -2 aβ vh S 0 v (N 0 h ) 2 (η h w 2 + w 3 ) w 1 v 2 2 aβ * hv N 0 h (η v w 6 + w 7 ) (w 2 + w 3 + w 4 ) h (η h w 2 + w 3 ) w 7 = ζ 1ζ 2 ,

  2 3 k 2 4 k 2 10 k 2 11 a 4 µ 2 h (N 0 v ) 2 β 2 vh , α 0 = k 2 3 k 2 4 k 2 10 a 2 µ 4 h (N 0 h ) 2 β 2 vh , and α 1 = 2k 3 k 4 k 10 k 11 a 2 µ 2 h N 0 h N 0 v β vh [2k 2 (k 10 aµ h β vhk 8 δγ h )k 3 k 4 k 10 aµ h β vh ] .Now we compute the discriminant ∆ := α 2 1 -4α 2 α 0 , to obtain:∆ = -16k 2 3 k 2 4 k 2 k 2 10 k 2 11 δγ h a 4 µ 4 h (N 0 h ) 2 (N 0 v ) 2 β 2 vh (k 10 aµ h β vh + k 8 k 2 ) (k 10 aµ h β vhδγ h k 8 ) .Equation (72) admits a real solution if and only if ∆ ≥ 0. This condition is equivalent to ψ := k 10 aµ h β vhδγ h k 8 ≤ 0 (73)

Acknowledgments

The first author (Hamadjam ABBOUBAKAR) thanks the Direction of UIT of Ngaoundere for the financial help granted under research missions in the year 2015. The first author also thanks the Department of Mathematics of ENS of Yaounde for their hospitality during the research visit where the work was initiated.