Mathematical Modelling and Optimal Control Technique and Analysis of Insecticides Susceptibility in Status on Mosquitoes and Community Knowledge on Malaria

1. Introduction

Malaria remains a significant public health challenge, particularly in tropical and subtropical regions where mosquitoes are prevalent. Effective management of this disease requires a comprehensive approach that combines mathematical modelling with optimal control strategies to enhance insecticide susceptibility in mosquito populations. In addition, raising community awareness of the dynamics of malaria transmission is crucial in reducing infection rates and promoting preventive measures. According to the World Malaria Report in [1], malaria continues to inflict unacceptable levels of illness and mortality. The most recent study estimates that as of April 18-20, 2023, there were 608,000 fatalities and 249 million cases worldwide. Since malaria is preventable and treated, lowering the burden of the disease and the death rate while maintaining the long-term goal of eliminating malaria should be a top priority for the entire world. Plasmodium parasites are the cause of malaria, and female Anopheles mosquitoes carry the disease. P. falciparum is the most harmful of the four species of malaria that affect humans. P. falciparum, P. vivax, P. malariae, and P. ovale. Of these, P. falciparum and P. vivax are the most common.Human infections with the zoonotic plasmodium P. knowlesi are also known to occur. Despite these advances, endemic malaria continues to exist in the six WHO areas, with the African Region bearing the brunt of the disease burden; an estimated 90% of all malaria deaths occur there. Two countries, Nigeria and the Democratic Republic of the Congo, account for about 40% of the global malarial fatality rate projected. Millions of people around the world lack access to malaria prevention and treatment, and the majority of cases and deaths go unreported. As the world population grows by 2030, more people will live in countries where malaria is a risk, further taxing the finances and health systems of the national malaria program [2].

Moreover, mosquito species like Anopheles and Culex are able to proliferate as they develop from larvae to adults due to the chemical insecticides (adulticides and lecticides) inadequate targeting of breeding habitats and residential areas [3,4,5]. Furthermore, parasites that feed on blood supply and guarantee hatchability are seen in female Anopheles gambiae and Culex quinquefasciatus mosquitoes. Furthermore, according to [6], disease is transmitted via vectors because of their increasing resistance to common pesticides. Additionally, studies on community host awareness and chemical pesticide sensitivity status have been conducted worldwide [7]. Despite the fact that no such study has been published, these frequent diseases spread by mosquitoes are becoming more prevalent worldwide. This is the first study to examine the most effective ways to control mosquito susceptibility to three distinct chemical insecticides: the usual technique (media campaign), propoxur (carbamate), permethrin (pyrethroids), and malathion (organophosphates). However, only a small number of nations have conducted prior tests on mosquitoes worldwide using the pyrethroid chemical insecticide components deltamethrin and permethrin. Nonetheless, the pesticides were suggested to be considered for the current investigation because of their resistance and long testing cycles. The dynamics of optimisation approaches are also used in this study to provide the model formulation in terms of optimal control theory. In this situation, the best control method provides powerful tools for creating and evaluating control plans in many scenarios. According to [8], there is evidence that malaria cases tend to cluster more when transmission levels are decreasing and get closer to zero. People who are exposed at the same time and place, like through a common vocation or shared travel to endemic areas, may cluster geographically, in small areas like households and neighbourhoods, or socially [9]. Malaria transmission at the community level may be decreased if clusters can be located and effectively targeted with interventions. Furthermore, in order to manage and measure mosquitoes and consequently reduce the transmission of related diseases, Hafez and Abbas, (2021) in [10] investigated insecticide resistance to insect growth regulators in Saudi Arabia. Over 700,000 people die each year from vector-borne infectious diseases, of which over 400,000 are from malaria alone. These infections are a major cause of death worldwide. Anopheles gambiae is a common mosquito species that is the main vector for human malaria transmission (Helmi, 2024) in [11]. In order to create a mathematical model, [12] combined three control strategies: mass awareness, treatment, and vector control. According to the results, the control of dengue fever is more effectively achieved when vector control, treatment, and public awareness are combined than when these interventions are used alone or in combination. Additionally, [13] examined the impact of abiotic parameters and species associations on the abundance of seven mosquito species in Spain’s Donana National Park. They evaluated the effects of abiotic parameters and species-to-species associations, which acted as stand-ins for species interactions, and developed three models with different parameters. A fractional order model for the transmission dynamics of malaria was constructed, according to [14], and it included two control strategies: the use of pesticides and health education campaigns. The findings indicated that the population’s exposure to malaria had been considerably reduced.

Thus, [15] developed a mathematical model to investigate the dynamics of malaria transmission in the presence of parasites and mosquitoes resistant to antimalarial drugs. A mathematical model for analysing the dynamics of malaria transmission is developed by [16]. In addition to affecting illness outcomes and healthcare systems, it takes into consideration consequences such as severe anaemia and organ failure. For malaria therapies to be effective, the results highlight the necessity of better vector management and complication control. A mathematical model of malaria transmission dynamics that is non-linear and deterministic was also proposed and examined by [17]. The ideal model’s results indicated that integrated control measures are superior to a single intervention for the eradication of malaria.[18] creates a compartmental model to assess how early and late treatment interventions affect the spread of malaria in children under $5years$ in order to inform efficient control measures. To examine the dynamics of malaria illness transmission and pesticide control strategies, [19] developed a mathematical model. The results were presented in a visual format. Insecticide spray has been shown to significantly affect the spread of malaria. Taking into consideration time-dependent treatment, immunisation, and ambient sanitation measures, Olaniyi, et al., (2025) in [20] developed a novel mathematical model for RVF. The results of the study shed light on the long-term dynamics of RVF in the community and provide effective preventative and control strategies with low intervention costs. To further reduce transmission between human populations and mosquitoes, Baroudi, em et al., (2025) in [21], provide an ideal approach that consists of health interventions, safety precautions, and awareness campaigns in dengue endemic areas. [22] created a model in mathematics. Numerical simulations were used to validate the analytical findings. According to the study’s findings, combining an ideal control plan with social media awareness efforts is the most economical way to treat malaria. Moreover, [23] developed a mathematical model of the best control method for biological and chemical management of mosquito populations. The results show two options with the best cost-effectiveness metrics and provide some useful information for their possible use in practical situations. A compartmental model was developed by [24] to demonstrate the spread of co-infection between cholera and malaria. The vector control method is the most efficient strategy to minimise the amount of new interactions, as the numerical results show. The application of a residual insecticide to the interior surfaces of walls, ceilings, windows, and doors is known as indoor residual spraying (IRS), and it is used to eliminate mosquitoes that are at rest and lessen the spread of malaria. Prior to the spread of malaria, the IRS, lasting nets insecticide (LNI), and traditional techniques are typically applied in the form of a campaign over a substantial region or a higher-risk region. Thus, it is advised that IRS, LNI, and traditional techniques be used in the homes of a confirmed case and their neighbours at roughly the same time [2].

The study is divided into a number of components, beginning with an introduction outlining the significance and goals of the paper. A thorough methodology section that describes the mathematical model and control techniques used to evaluate pesticide susceptibility comes next. Findings from the modelling analysis and community surveys are presented in the results section, and their interpretation in relation to malaria awareness and control is discussed.

The aim of the paper is to develop a mathematical model that predicts how insecticide susceptibility in mosquito populations affects malaria transmission dynamics. Additionally, the study seeks to assess the impact of community awareness initiatives on reducing the spread of malaria. Through this approach, the paper intends to provide insights for creating more effective vector control strategies and public health interventions.

2. Materials and Methods

The model is formulated by incorporating variables that represent the insecticide susceptibility status of mosquito populations and the level of community awareness about malaria prevention. Differential equations are used to describe the dynamics of mosquito populations, taking into account factors such as birth and death rates, as well as the impact of insecticide use. Additionally, parameters related to community education initiatives are included to assess how increased awareness can influence the effectiveness of malaria control strategies.

Therefore, the diagram can also be used to identify areas where controls may be most effective, which can be seen in Figure 1 as follow;

Figure 1. Diagram of a Malaria Transmission. Source: Authors.

Figure 1. Diagram of a Malaria Transmission. Source: Authors.

Table 1. Description of Variables.

Variables	Biological description
$S_h\left(t\right)$	The number susceptible humans with time
$A_h\left(t\right)$	The number of aware humans with time
$E_h\left(t\right)$	The number of exposed humans with time
$I_h\left(t\right)$	The infectious humans with time
$R_h\left(t\right)$	The recovered humans with time
$S_m\left(t\right)$	The number of adult healthy mosquitoes with time
$E_m\left(t\right)$	The number of exposed mosquitoes with time
$I_m\left(t\right)$	The number of infectious mosquitoes with time
$A_q\left(t\right)$	The number of Aquatic mosquitoes (eggs, larvae, pupae) with time

Table 1. Description of Variables.

Variables	Biological description
$S_h\left(t\right)$	The number susceptible humans with time
$A_h\left(t\right)$	The number of aware humans with time
$E_h\left(t\right)$	The number of exposed humans with time
$I_h\left(t\right)$	The infectious humans with time
$R_h\left(t\right)$	The recovered humans with time
$S_m\left(t\right)$	The number of adult healthy mosquitoes with time
$E_m\left(t\right)$	The number of exposed mosquitoes with time
$I_m\left(t\right)$	The number of infectious mosquitoes with time
$A_q\left(t\right)$	The number of Aquatic mosquitoes (eggs, larvae, pupae) with time

Table 2. Description of Parameters.

Parameters	Biological description
${\Lambda }_h$	Recruitment per-capita rate of susceptible humans
${\Lambda }_m$	Recruitment per-capital rate of susceptible mosquitoes
$\epsilon$	Rate of aware human return to susceptible humans
$\alpha$	Rate of recovered human become aware humans
${\mu }_h$	Natural death rates of all humans
${\mu }_m$	Natural death rates of all mosquitoes
${\rho }_h$	Transmission rates per-capita of humans
${\rho }_m$	Transmission rates per-capita of mosquitoes
${\beta }_h$	Per- capita contact rates of human with infected mosquitoes
${\beta }_m$	Per- capita contact rates of healthy mosquitoes with infectious humans
a	Rates of exposed human move to infectious humans
b	Rate of recovered humans
$\theta$	Progression rate of exposed mosquitoes become infectious mosquitoes
${\delta }_h$	Disease induced-mortality rate of infectious human
${\rho }_q$	Proportion rate of the oviposition
${\xi }_m$	Proportion rate of non-infected eggs laid by infected mosquitoes.
$\delta$	Proportion rate in which mosquitoes mature
${\delta }_m$	Mortality rate of infectious mosquitoes due to human activities
${\delta }_a$	Mortality rate of aquatic mosquitoes due to human activities

Table 2. Description of Parameters.

Parameters	Biological description
${\Lambda }_h$	Recruitment per-capita rate of susceptible humans
${\Lambda }_m$	Recruitment per-capital rate of susceptible mosquitoes
$\epsilon$	Rate of aware human return to susceptible humans
$\alpha$	Rate of recovered human become aware humans
${\mu }_h$	Natural death rates of all humans
${\mu }_m$	Natural death rates of all mosquitoes
${\rho }_h$	Transmission rates per-capita of humans
${\rho }_m$	Transmission rates per-capita of mosquitoes
${\beta }_h$	Per- capita contact rates of human with infected mosquitoes
${\beta }_m$	Per- capita contact rates of healthy mosquitoes with infectious humans
a	Rates of exposed human move to infectious humans
b	Rate of recovered humans
$\theta$	Progression rate of exposed mosquitoes become infectious mosquitoes
${\delta }_h$	Disease induced-mortality rate of infectious human
${\rho }_q$	Proportion rate of the oviposition
${\xi }_m$	Proportion rate of non-infected eggs laid by infected mosquitoes.
$\delta$	Proportion rate in which mosquitoes mature
${\delta }_m$	Mortality rate of infectious mosquitoes due to human activities
${\delta }_a$	Mortality rate of aquatic mosquitoes due to human activities

2.1. Model Assumption

The model assumes that the mosquito population is homogeneously mixed, meaning that each mosquito has an equal chance of coming into contact with an insecticides. It also presumes that community awareness efforts are uniformly distributed and have a consistent impact on reducing malaria transmission. Furthermore, the model considers environmental factors to be constant over time, allowing for a focus on the interactions between control measures and mosquito behaviour. It assumed that the human and mosquito populations are divided into separate classes, with each distinct one assigned to only one class at a time. In this study, mosquitoes go through an aquatic (water-based) stage before becoming adults with the rate ${\displaystyle \frac{1}{\delta }}$. Only adult mosquitoes are capable of transmitting malaria. It is assumed that malaria spreads through bites?susceptible humans can get infected by bites from infectious mosquitoes, and mosquitoes can become infectious after biting infected humans. Also assumed that the people who recover from malaria may become more aware of the disease and help increase overall community awareness. Therefore, insecticide use increases the death rate of mosquitoes, especially those that are infected. Finally, the community awareness helps reduce the likelihood of people getting infected by promoting protective actions and behaviours.

2.2. Model Description

Susceptible humans $S_h\left(t\right)$; are humans who have not yet been infected with malaria but are at risk of contracting the disease. These humans are increasing the population through immigration and birth rate. Then decreasing the population through the infectious mosquito rate and the natural death rate of mosquitoes, which gives;

$${\displaystyle \frac{d{S}_h}{dt}}={\Lambda }_h+\epsilon A_h-{\rho }_h{\beta }_h{S}_h{I}_m-{\mu }_h{S}_h.$$

Aware humans $A_h\left(t\right)$;these population increasing through recovered rate of human and decreasing to susceptible human class. This class are knowledgeable about the life cycle of mosquitoes and the role they play in malaria transmission, thus;

$${\displaystyle \frac{d{A}_h}{dt}}=\alpha R_h-(\epsilon +{\mu }_h)A_h.$$

Exposed humans $E_h\left(t\right)$; are individuals who have been bitten by mosquitoes carrying the malaria parasite but have not yet developed symptoms of the disease. This class increasing through contact rates of human with mosquitoes and decreasing with some rates move to infectious human, gives;

$${\displaystyle \frac{d{E}_h}{dt}}={\rho }_h{\beta }_h{S}_h{I}_m-({\mu }_h+a)E_h.$$

Infectious humans $I_h\left(h\right)$; play a crucial role in the transmission dynamics of malaria, as they serve as hosts for the malaria parasites. When a mosquito bites an infectious human, it can acquire the parasites, which then develop within the mosquito before being transmitted to other humans, thus;

$${\displaystyle \frac{d{I}_h}{dt}}=a{E}_h-({\mu }_h+b+{\delta }_h)I_h.$$

Recovered human $R_h\left(t\right)$; a recovered human in this context refers to an individual who has successfully overcome a malaria infection and is now immune to the specific strain they encountered. thus;

$${\displaystyle \frac{d{R}_h}{dt}}=b{I}_h-({\mu }_h+\alpha )R_h.$$

Susceptible mosquitoes $S_m\left(t\right)$; are those that are vulnerable to the effects of insecticides, meaning they can be effectively controlled or eliminated through chemical interventions, gives;

$${\displaystyle \frac{d{S}_m}{dt}}={\Lambda }_m+\delta A_g-{\rho }_m{\beta }_m{S}_m{I}_h-{\mu }_m{S}_m.$$

Exposed mosquitoes$E_m\left(t\right)$; are those that have come into contact with insecticides but have not yet succumbed to their effects, thus;

$${\displaystyle \frac{d{E}_m}{dt}}={\rho }_m{\beta }_m{S}_m{I}_h-({\mu }_m+\theta )E_m.$$

Infectious mosquitoes $I_m\left(t\right)$ are produced by the maturation of infected aquatic mosquitoes and the infection of susceptible mosquitoes. thus;

$${\displaystyle \frac{d{I}_m}{dt}}=\theta E_m-({\mu }_m+{\delta }_m)I_m.$$

Aquatic mosquitoes $A_q\left(t\right)$; increasing via oviposition by $S_m\left(t\right)$ and infectious $I_m\left(t\right)$ mosquitoes at ${\rho }_q$ and ${\xi }_q{\rho }_q$ respectively, where ${\rho }_q$ is the oviposition rate while ${\xi }_q$ is a proportion of non-infected eggs laid by infected mosquitoes. It decreases through mature mosquitoes at a rate $\delta$, die naturally at ${\mu }_m$ and mortality rate at $\delta$, thus;

$${\displaystyle \frac{d{A}_q}{dt}}={\rho }_q{\beta }_m(S_m+{\xi }_q{I}_m)-\delta A_q-({\mu }_m+{\delta }_a)A_q.$$

The following non-linear differential equations represent the model description for the transmission of mosquito diseases in the community and correspond to the model diagram in Figure 1, gives;

$$\begin{array}{c}\hfill \left.\begin{array}{c}{\displaystyle \frac{d{S}_h}{dt}}={\Lambda }_h+\epsilon A_h-{\rho }_h{\beta }_h{S}_h{I}_m-{\mu }_h{S}_h,\hfill \\ \\ {\displaystyle \frac{d{A}_h}{dt}}=\alpha R_h-(\epsilon +{\mu }_h)A_h,\hfill \\ \\ {\displaystyle \frac{d{E}_h}{dt}}={\rho }_h{\beta }_h{S}_h{I}_m-({\mu }_h+a)E_h,\hfill \\ \\ {\displaystyle \frac{d{I}_h}{dt}}=a{E}_h-({\mu }_h+b+{\delta }_h)I_h,\hfill \\ \\ {\displaystyle \frac{d{R}_h}{dt}}=b{I}_h-({\mu }_h+\alpha )R_h,\hfill \\ \\ {\displaystyle \frac{d{S}_m}{dt}}={\Lambda }_m+\delta A_q-{\rho }_m{\beta }_m{S}_m{I}_h-{\mu }_m{S}_m,\hfill \\ \\ {\displaystyle \frac{d{E}_m}{dt}}={\rho }_m{\beta }_m{S}_m{I}_h-({\mu }_m+\theta )E_m,\hfill \\ \\ {\displaystyle \frac{d{I}_m}{dt}}=\theta E_m-({\mu }_m+{\delta }_m)I_m,\hfill \\ \\ {\displaystyle \frac{d{A}_q}{dt}}={\rho }_q{\beta }_m(S_m+{\xi }_m{I}_m)-\delta A_q-({\mu }_m+{\delta }_a)A_q.\hfill \end{array}\right\}\end{array}$$

(1)

with the initial conditions;

$S_h\left(0\right)\ge 0,A_h\left(0\right)\ge 0,E_h\left(0\right)\ge 0,I_h\left(0\right)\ge 0,R_h\left(0\right)\ge 0,S_m\left(0\right)\ge 0,E_m\left(0\right)\ge 0,I_m\left(0\right)\ge 0$, and $A\left(q\right)\left(0\right)\ge 0$.

2.3. Basic Properties of the Model

The basic properties of the model include the parameters that govern the interaction between mosquitoes and humans, such as transmission rates and insecticide effectiveness. The feasible region defines the set of all possible states of the system, ensuring that population sizes remain non-negative and within biologically realistic limits. The invariant region refers to a subset of the feasible region where the system’s dynamics are constrained to remain over time, reflecting stable population levels and effective malaria control.

2.3.1. Feasible and Invariant Region

The feasible region in this context refers to the set of all possible solutions that satisfy the constraints of the model. Therefore, using the idea of Haile, et al., (2025) in [25] as system (1) will be examined and divided into two, which are both human and mosquito classes, respectively, and given as $\Omega ={\Omega }_{h}X{\Omega }_m.$

Theorem 1.

Let system (1) be the set of $\left\{S_h,A_h,E_h,I_h,R_h,S_m,E_m,I_m,A_q\right\}\in R_{+}^9,$ then the feasible region Ω contain the all solution.

Proof of Theorem 1.

Suppose $t\ge 0$, then $\left\{S_h,A_h,E_h,I_h,R_h,S_m,E_m,I_m,A_q\right\}\in R_{+}^9$. To examined the dynamics of system (1), the $\Omega$ is positively invariant. This region is determined by considering the entire human population. $S_h\left(t\right)+A_h\left(t\right)+E_h\left(t\right)+I_h\left(t\right)+R_h\left(t\right)$ = $N_h\left(t\right)$ Combining the system (1) of the first five equations and differentiating both sides with respect to t yields;

$$\begin{array}{c}\hfill {\displaystyle \frac{d}{dt}}\left(S_h+A_h+E_h+I_h+R_h\right)={\Lambda }_h-{\mu }_h{N}_h-\delta I_h.\end{array}$$

(2)

This implies that;

${\displaystyle \frac{d}{dt}}\left(S_h+A_h+E_h+I_h+R_h\right)\le {\Lambda }_h-{\mu }_h{N}_h.$ Thus;

$${\displaystyle \frac{d{N}_h}{dt}}\le {\Lambda }_h-{\mu }_h{N}_h.$$

(3)

Integrating and simplifying both side of equation (3), obtain;

$${\Lambda }_h-{\mu }_h{N}_h\ge B_{1}exp(-\mu t).$$

(4)

Where $B_1$ is constant value. Using the initial condition, rearrange equation (4), give;

$$N_h\le {\displaystyle \frac{{\Lambda }_h}{{\mu }_h}}-{\displaystyle \frac{N_{h0}exp(-\mu t)}{{\mu }_h}}.$$

(5)

Thus, as $t⟶\infty$, the total human population $N_h⟶{\Lambda }_h$, which indicates that $0\le N_h\le {\displaystyle \frac{{\Lambda }_h}{{\mu }_h}}.$ The invariant region of system (1) for the human population therefore yields:

$${\Omega }_h=\left\{\left(S_h,A_h,E_h,I_h,R_h\right)\in R_{+}^5:0\le S_h+A_h+E_h+I_h+R_h\le {\displaystyle \frac{{\Lambda }_h}{{\mu }_h}}.\right\}$$

(6)

Thus, is positively invariant. Subsequently, the total population of mosquito in system (1) follow as;

$$N_m\left(t\right)=S_m+E_m+I_m+A_m.$$

(7)

Also, differentiating both side of equation (7) w.r.t, t obtain;

$${\displaystyle \frac{d}{dt}}\left(S_m+E_m+I_m+A_q\right)={\Lambda }_m-{\mu }_m{N}_h-{\delta }_m{I}_m-{\delta }_a{A}_q.$$

(8)

$${\displaystyle \frac{d{N}_m}{dt}}\le {\Lambda }_m-{\mu }_m{N}_m.$$

(9)

Solving equation (9), gives $0\le N_m\le {\displaystyle \frac{{\Lambda }_m}{{\mu }_m}}.$ Hence, the invariant region of system (1) for mosquito population gives;

$${\Omega }_m=\left\{\left(S_m,E_m,I_m,A_q\right)\in R_{+}^4:0\le S_m+E_m+A_q\le {\displaystyle \frac{{\Lambda }_m}{{\mu }_m}}.\right\}$$

(10)

Therefore, the invariant region of the whole system (1) is given as;

$$\Omega ={\Omega }_{h}X{\Omega }_m=\left\{\left(S_h,A_h,E_h,I_h,R_h,S_m,E_m,I_m,A_q\right)\in R_{+}^9:N_h\le {\displaystyle \frac{{\Lambda }_h}{{\mu }_h}},N_m\le {\displaystyle \frac{{\Lambda }_m}{{\mu }_m}}.\right\}$$

(11)

The prove is complete and is positively invariant. Finally, all the solution set of system (1) is bounded in $\Omega .$ □

2.3.2. Positivity of the Model Solutions

Suppose $t\ge 0$, for all solution to system (1) with positive initial conditions will stay positive .

Theorem 2.

System (1) can be solved as follows, given positive initial conditions; $S_h,A_h,E_h,I_h,R_h,S_m,E_m,I_m$ and $A_q$ with positive initial conditions. Also, $S_h\left(0\right),A_h\left(0\right),E_h\left(0\right),I_h\left(0\right),R_h\left(0\right),S_m\left(0\right),E_{m\left(0\right)},I_m\left(0\right)$ and $A_q\left(0\right)$ remains positive for all $t\ge 0$.

Proof of Theorem 2.

Assume that every state variable in system (1) is positive. Starting with (1), the first equation of system (1) is ${\displaystyle \frac{d{S}_h}{dt}}={\Lambda }_h+\epsilon A_h-{\rho }_h{\beta }_h{S}_h{I}_m-{\mu }_h{S}_h$,

this implies;

$${\displaystyle \frac{d{S}_h}{dt}}\ge -\left({\rho }_h{\beta }_h{I}_m-{\mu }_h\right)S_h.$$

(12)

Integrating by applying the method of separation of variable along with initial condition, gives;

$$S_h\left(t\right)\ge S_h\left(0\right)exp-\left({\rho }_h{\beta }_h{I}_m-{\mu }_h\right)t\ge 0.$$

(13)

Also, in the second equation of system (1) as ${\displaystyle \frac{d{A}_h}{dt}}=\alpha R_h-(\epsilon +{\mu }_h)A_h$; this implies;

$${\displaystyle \frac{d{A}_h}{dt}}\ge -(\epsilon +{\mu }_h)A_h$$

(14)

Integrating the both side and applying initial condition, gives;

$$A_h\left(t\right)\ge A_h\left(0\right)exp-\left(\epsilon +{\mu }_h\right)t\ge 0.$$

(15)

The third equation in the system (1) , as ${\displaystyle \frac{d{E}_h}{dt}}={\rho }_h{\beta }_h{S}_h{I}_m-({\mu }_h+a)E_h$, and implies;

$${\displaystyle \frac{d{E}_h}{dt}}\ge -({\mu }_h+a)E_h$$

(16)

Integrating the both side and applying initial condition, gives;

$$E_h\left(t\right)\ge E_h\left(0\right)exp-\left(a+{\mu }_h\right)t\ge 0.$$

(17)

Furthermore, by applying initial condition on all the equations of system (1), obtains;

$$I_h\left(t\right)\ge I_h\left(0\right)exp-\left(b+{\mu }_h+{\delta }_h\right)t\ge 0.$$

(18)

$$R_h\left(t\right)\ge R_h\left(0\right)exp-\left(\alpha +{\mu }_h\right)t\ge 0.$$

(19)

$$S_m\left(t\right)\ge S_m\left(0\right)exp-\left({\rho }_m{\beta }_m{I}_h+{\mu }_m\right)t\ge 0.$$

(20)

$$E_m\left(t\right)\ge E_m\left(0\right)exp-\left({\mu }_m+\theta \right)t\ge 0.$$

(21)

$$I_m\left(t\right)\ge I_m\left(0\right)exp-\left({\mu }_m+{\delta }_m\right)t\ge 0.$$

(22)

$$A_q\left(t\right)\ge A_q\left(0\right)exp-\left(\delta +{\mu }_m+{\delta }_a\right)t\ge 0.$$

(23)

Thus, the system (1) is mathematically well posed (Hethcote, 2000) in [26]. This complete the prove. □

2.4. Mosquito Disease Free Equilibrium Points

The Mosquito Disease Free Equilibrium Point (DFEP) refers to a state where the mosquito population is present, but there are no cases of the disease being transmitted within the community, meaning $E_h=I_h=E_m=I_m=A_q=0$. The disease free equilibrium in the system (1) is determined using Maple23 software to identify and set all the model equations to zero.Thus, the disease-free equilibrium point of system (1) is given as $E_{*}=\left(S_{h*},A_{h*},E_{h*},I_{h*},R_{h*},S_{m*},E_{m*},I_{m*},A_{q*}\right)$. Also, this implies that the disease-free equilibrium is given by;

$$E_{*}=\left({\displaystyle \frac{{\Lambda }_h}{{\mu }_h}},0,0,0,0,{\displaystyle \frac{{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right)}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}},0,0,{\displaystyle \frac{{\Lambda }_m{\rho }_q{\beta }_m}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}\right).$$

(24)

2.5. Endemic Equilibrium Point

At the endemic equilibrium, the population dynamics of mosquitoes and the incidence of malaria cases reach a steady state. This balance reflects a consistent rate of transmission and infection, where the number of new cases is equal to the number of recoveries or deaths. Modelling this equilibrium is crucial for designing effective control strategies that can reduce the prevalence of malaria in affected communities and is denoted as $E^{*}$.

Therefore, EEP are as follows;

$$\begin{array}{c}{S}_h^{*}={\displaystyle \frac{A_7{A}_8\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)\left(ab\alpha A_1+A_2{A}_3{A}_4{A}_5\right)}{{\beta }_m{\mu }_h{A}_2{A}_5{\rho }_{m}a\left(\theta A_6{\rho }_q{\xi }_m{\beta }_m+A_7{A}_8{A}_9\right)}},\\ \\ A_h^{*}={\displaystyle \frac{b\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)A_7{A}_8\alpha }{\left(\theta A_6{\rho }_q{\xi }_m{\beta }_m+A_7{A}_8{A}_9\right){\rho }_m{\beta }_m{A}_2{A}_5}},\\ \\ E_h^{*}={\displaystyle \frac{A_4\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)A_7{A}_8}{\left(\theta A_6{\rho }_q{\xi }_m{\beta }_m+A_7{A}_8{A}_9\right){\rho }_m{\beta }_{m}a}},\\ \\ I_h^{*}={\displaystyle \frac{\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)A_7{A}_8}{\left(\theta A_6{\rho }_q{\xi }_m{\beta }_m+A_7{A}_8{A}_9\right){\rho }_m{\beta }_m}},\\ \\ R_h^{*}={\displaystyle \frac{b\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)A_7{A}_8}{\left(\theta A_6{\rho }_q{\xi }_m{\beta }_m+A_7{A}_8{A}_9\right){\rho }_m{\beta }_m{A}_5}},\\ \\ S_m^{*}={\displaystyle \frac{A_3{A}_4{\mu }_h{A}_2{A}_5\left(\theta A_6{\rho }_q{\xi }_m{\beta }_m+A_7{A}_8{A}_9\right)}{\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)\left(ab\alpha A_1+A_2{A}_3{A}_4{A}_5\right){\rho }_h{\beta }_h\theta }},\\ \\ E_m^{*}={\displaystyle \frac{A_8{\mu }_h{A}_2{A}_5{A}_3{A}_4}{{\rho }_h{\beta }_h\left(ab\alpha A_1-A_2{A}_3{A}_4{A}_5\right)\theta }},\\ \\ I_m^{*}={\displaystyle \frac{A_3{A}_4{\mu }_h{A}_2{A}_5}{{\rho }_h{\beta }_h\left(ab\alpha A_1-A_2{A}_3{A}_4{A}_5\right)}},\\ \\ A_q^{*}={\displaystyle \frac{A_5{A}_2{\mu }_h{A}_4{A}_3{\rho }_q{\beta }_m\left(\theta {\mu }_m{\xi }_m+A_7{A}_8\right)}{\left(A_6{\rho }_q{\beta }_m+A_9{\mu }_m\right)\left(ab\alpha A_1+A_2{A}_3{A}_4{A}_5\right){\rho }_h{\beta }_h\theta }}.\end{array}$$

(25)

2.6. Reproduction Number

The Basic Reproduction Number, often denoted as $\left(R_0\right)$, is a key epidemiological metric used to describe the contagiousness or transmissibility of infectious agents. It represents the average number of secondary cases generated by one primary case in a completely susceptible population. Understanding $\left(R_0\right)$ helps in assessing the potential spread of malaria and the effectiveness of control measures. The following result is produced by standard method in [27,28].

Considering the infected classes, which are $E_h,I_h,E_m$ and $I_m$, gives;

$f=\left[\begin{array}{c}{\rho }_h{\beta }_h{S}_h\\ 0\\ {\rho }_m{\beta }_m{S}_m\\ 0\end{array}\right]$ and $v=\left[\begin{array}{c}(a+{\mu }_h{E}_h)\\ -a{E}_h+(b+{\mu }_h+{\delta }_h{I}_h)\\ ({\mu }_m+\theta )E_m\\ -\theta E_m+({\mu }_m+{\delta }_m{I}_m)\end{array}\right]$.

Solving the Jacobian of matrices f and v, then differentiating w.r.t $E_h,I_h,E_m$ and $I_m$, gives;

$F=\left[\begin{array}{cccc}0& 0& 0& {\rho }_h{\beta }_h{\displaystyle \frac{{\Lambda }_h}{{\mu }_h}}\\ 0& 0& 0& 0\\ 0& {\rho }_m{\beta }_m{S}_{m*}& 0& 0\\ 0& 0& 0& 0\end{array}\right]$ and

$V=\left[\begin{array}{cccc}(a+{\mu }_h)& 0& 0& 0\\ -a& (b+{\mu }_h+{\delta }_h)& 0& 0\\ 0& 0& ({\mu }_m+\theta )& 0\\ 0& 0& -\theta & ({\mu }_m+{\delta }_m)\end{array}\right]$.

Therefore, using the method in [29,30], the result follows as;

$V^{-1}=\left[\begin{array}{cccc}{\displaystyle \frac{1}{a+{\mu }_h}}& 0& 0& 0\\ {\displaystyle \frac{a}{\left(b+{\mu }_h+{\delta }_h\right)\left(a+{\mu }_h\right)}}& {\displaystyle \frac{1}{b+{\mu }_h+{\delta }_h}}& 0& 0\\ 0& 0& {\displaystyle \frac{1}{{\mu }_m+\theta }}& 0\\ 0& 0& {\displaystyle \frac{\theta }{\left({\mu }_m+\theta \right)\left({\mu }_m+{\delta }_m\right)}}& {\displaystyle \frac{1}{{\mu }_m+{\delta }_m}}\end{array}\right].$

Subsequently;

$F{V}^{-1}=\left[\begin{array}{cccc}0& 0& {\displaystyle \frac{{\rho }_h{\beta }_h{\Lambda }_h\theta }{{\mu }_h\left({\mu }_m+\theta \right)\left({\mu }_m+{\delta }_m\right)}}& {\displaystyle \frac{{\rho }_h{\beta }_h{\Lambda }_h}{{\mu }_h\left({\mu }_m+{\delta }_m\right)}}\\ 0& 0& 0& 0\\ {\displaystyle \frac{{\rho }_m{\beta }_m{S}_{m*}a}{\left(b+{\mu }_h+{\delta }_h\right)\left(a+{\mu }_h\right)}}& {\displaystyle \frac{{\rho }_m{\beta }_m{S}_{m*}}{b+{\mu }_h+{\delta }_h}}& 0& 0\\ 0& 0& 0& 0\end{array}\right].$

Therefore, the basic reproduction number of the model system (1) is denoted by $R_0$ = $\rho F{V}^{-1}$, where $\rho$ is the spectral radius of the product $F{V}^{-1}$, which gives;

$$R_0=\sqrt{{\displaystyle \frac{{\mu }_h\left(a+{\mu }_h\right)\left(b+{\mu }_h+{\delta }_h\right)\left({\mu }_m+\theta \right)\left({\mu }_m+{\delta }_m\right){\rho }_m{\beta }_m{S}_{m*}a{\rho }_h{\beta }_h{\Lambda }_h\theta }{{\mu }_h\left(a+{\mu }_h\right)\left(b+{\mu }_h+{\delta }_h\right)\left({\mu }_m+\theta \right)\left({\mu }_m+{\delta }_m\right)}}}.$$

(26)

Where;

$S_{m*}={\displaystyle \frac{{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right)}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}$.

Therefore, by further simplification, the basic reproduction number $R_0$, gives;

$$R_0=\sqrt{{\displaystyle \frac{{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right){\rho }_m{\beta }_{m}a{\rho }_h{\beta }_h{\Lambda }_h\theta }{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}}.$$

(27)

2.7. Local Stability of Disease Free Equilibrium

The local stability of the disease-free equilibrium is crucial in understanding how effective control measures can be in eradicating malaria. A disease-free equilibrium is stable when the basic reproduction number, $R_0$, is less than one, indicating that the infection will eventually die out. Mathematical models can help identify conditions under which this stability is achieved, guiding strategies for insecticide use and community education efforts. The following theorem analyses the local stability of DFE.

Theorem 3.

The DFE of system (1), denoted by $E_{*}$, is unstable if $R_0>1$ and locally asymptotically stable (LAS) in Ω if $R_0<1$.

Proof of Theorem 3.

At disease-free equilibrium $E_{*}$ the Jacobian matrix $J\left(E_{*}\right)$ is gives as;

$$J\left(E_{*}\right)=\left[\begin{array}{ccccccccc}-{\mu }_h& \epsilon & 0& 0& 0& 0& 0& -B& 0\\ 0& A_1& 0& 0& \alpha & 0& 0& 0& 0\\ 0& 0& A_{2}a& 0& 0& 0& 0& B& 0\\ 0& 0& a& A_3-& 0& 0& 0& 0& 0\\ 0& 0& 0& b& A_4& 0& 0& 0& 0\\ 0& 0& 0& -A_8& 0& {\mu }_m& 0& 0& \delta \\ 0& 0& 0& A_8& 0& 0& A_5& 0& 0\\ 0& 0& 0& 0& 0& 0& \theta & A_6& 0\\ 0& 0& 0& 0& 0& {\rho }_q{\beta }_m& 0& A_9& A_7\end{array}\right].$$

(28)

The eigenvalues of the Jacobian matrix are given by finding the characteristic polynomial, which follow as; $\mid J\left(E_{*}\right)-\lambda I\mid =0$, obtain;

$$\left[\begin{array}{ccccccccc}-{\mu }_h-\lambda & \epsilon & 0& 0& 0& 0& 0& -B& 0\\ 0& A_1-\lambda & 0& 0& \alpha & 0& 0& 0& 0\\ 0& 0& A_2-\lambda & 0& 0& 0& 0& B& 0\\ 0& 0& a& A_3-\lambda & 0& 0& 0& 0& 0\\ 0& 0& 0& b& A_4-\lambda & 0& 0& 0& 0\\ 0& 0& 0& -A_8& 0& {\mu }_m-\lambda & 0& 0& \delta \\ 0& 0& 0& A_8& 0& 0& A_5-\lambda & 0& 0\\ 0& 0& 0& 0& 0& 0& \theta & A_6-\lambda & 0\\ 0& 0& 0& 0& 0& {\rho }_q{\beta }_m& 0& A_9& A_7-\lambda \end{array}\right]=0.$$

(29)

The focus is on $9\times 9$ matrix, by considering the rows and columns due to the presence of zeros in the reduction process of equation (29), gives;

$$\left[\begin{array}{ccccc}{A}_2-\lambda & 0& 0& B& 0\\ a& A_3-\lambda & 0& 0& 0\\ 0& -A_8& {\mu }_m-\lambda & 0& \delta \\ 0& 0& 0& A_6-\lambda & 0\\ 0& 0& {\rho }_q{\beta }_m& A_9& A_7-\lambda \end{array}\right]=0.$$

(30)

The characteristic polynomial of Jacobian matrix in equation (30), is given by;

$$\left(\lambda +A_6\right)\left(A_2+\lambda \right)\left(A_3+\lambda \right)\left({\lambda }^2+\left(A_7-{\mu }_m\right)\lambda -\delta {\rho }_q{\beta }_m+A_7{\mu }_m\right)=0.$$

(31)

Where;

$$\begin{gathered} A_1=-(\epsilon +{\mu }_h), \\ {A}_2=-(a+{\mu }_h), \\ {A}_3=-(b+{\mu }_h+{\delta }_h), \\ {A}_4=-(\alpha +{\mu }_h), \\ {A}_5=-({\mu }_m+\theta ), \\ {A}_6=-({\mu }_m+{\delta }_m), \\ {A}_7=-(\delta +{\mu }_m+{\delta }_a), \\ {A}_8={\displaystyle \frac{{\rho }_m{\beta }_m{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right)}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}, \\ {A}_9={\rho }_q{\beta }_m{\xi }_m\mathrm{and} \\ B={\displaystyle \frac{{\rho }_h{\beta }_h{\Lambda }_h}{{\mu }_h}}. \end{gathered}$$

Thus, from first expression of equation (31);

$\lambda =-({\mu }_m+{\delta }_m);{\lambda }_1=-(a+{\mu }_h)$ and ${\lambda }_2=-(b+{\mu }_h+{\delta }_h).$

Subsequently, the second expression of equation (31) can also be in the form of quadratic expression, which is;

$${\lambda }^2+N_1\lambda +N_2=0.$$

(32)

Where;

$N_1=-(\delta +{\mu }_m+{\delta }_a+{\mu }_m),$

$N_2=\delta {\rho }_q{\beta }_m-{\mu }_m(\delta +{\mu }_m+{\delta }_a).$

Applying the Routh-Harwitz criteria in [25,31] and by the principle of equation (31), all the real root of the equation has strictly negative and iff $N_1\ge 0,N_2\ge 0$, Then $N_1,N_2\ge 0$.This shown that $N_2\ge 0$ because it is the sum of positive parameters [25,31]. Hence, the DFE is locally asymptotically stable if $R_0<0$.

However, when $R_0\ge 0$, this shows that $N_1<0$ and by Descartes’s rule of signs in [25] there is exactly one sign change in the coefficient of the characteristic polynomial in equation (31). So, there is one eigenvalue with positive real part and DFE is unstable. Hence the proved. □

2.8. Global Stability of Disease Free Equilibrium

The global stability of the disease-free equilibrium is crucial for ensuring that malaria transmission can be effectively eradicated. The disease-free state remains stable despite potential perturbations. This involves ensuring the basic reproduction number is less than one, indicating that malaria cannot sustain itself within the community.

Theorem 4.

If $R_0<1$, the disease-free equilibrium for system (1) is globally asymptotically stable in the feasible region.

Proof of Theorem 4.

Using the idea of Haile, et al., (2025) in [25] to proof the theorem . The Lyapunov function method is then used to generate an appropriate Lyapunov function $V\left(t\right)$ following the procedure described in [25]. This will show the global asymptotic stability of the equilibrium point.

Let the Lyapunov function be $V=x_1{E}_h+x_2{I}_h+x_3{E}_m+x_4{I}_m$ , where $x_1,x_2,x_3$ and $x_4$ are non-negative constant and $E_h,I_h,E_m$ and $I_m$ non-negative infected classes. Then differentiate the Lyapunov w.r.t time, which gives;

$${\displaystyle \frac{dV}{dt}}=x_1{\displaystyle \frac{d{E}_h}{dt}}+x_2{\displaystyle \frac{d{I}_h}{dt}}+x_3{\displaystyle \frac{d{E}_m}{dt}}+x_4{\displaystyle \frac{I_m}{dt}}.$$

(33)

Putting the values of ${\displaystyle \frac{d{E}_h}{dt}},{\displaystyle \frac{d{I}_h}{dt}},{\displaystyle \frac{d{E}_m}{dt}}$ and ${\displaystyle \frac{I_m}{dt}}$ of system (1) to equation (33) by solving it and collecting like terms of the equation, obtain;

$$\begin{gathered} {\displaystyle \frac{dV}{dt}}=x_1\left[{\rho }_h{\beta }_h{S}_h{I}_m-(a+{\mu }_h)E_h\right]+x_2\left[a{E}_h-(b+{\mu }_h+{\delta }_h)I_h\right] \\ +x_3\left[{\rho }_m{\beta }_m{S}_m{I}_h-(\theta +{\mu }_m)E_m\right]+x_4\left[\theta E_m-({\mu }_m+{\delta }_m)I_m\right] \end{gathered}$$.

This implies;

$$\begin{gathered} {\displaystyle \frac{dV}{dt}}=x_1{\rho }_h{\beta }_h{S}_h{I}_m-x_1(a+{\mu }_h)E_h+x_{2}a{E}_h-x_2(b+{\mu }_h+{\delta }_h)I_h+x_3{\rho }_m{\beta }_m{S}_m{I}_h \\ -x_3(\theta +{\mu }_m)E_m+x_4\theta E_m-x_4({\mu }_m+{\delta }_m)I_m \end{gathered}$$.

Subsequently;

$$\begin{gathered} {\displaystyle \frac{dV}{dt}}=\left(x_3{\rho }_m{\beta }_m{S}_m-x_2(b+{\mu }_h+{\delta }_h)\right)I_h+\left(x_{2}a-x_1(a+{\mu }_h)\right)E_h \\ +\left(x_1{\rho }_h{\beta }_h{S}_h-x_4({\mu }_m+{\delta }_m)\right)I_m+\left(x_4\theta -({\mu }_m+\theta )\right)E_m. \end{gathered}$$

Taking the coefficient of $I_h,E_h$ and $I_m=0$ and obtain;

$${\displaystyle \frac{dV}{dt}}=\left(x_4\theta -x_3({\mu }_m+\theta )\right)E_m.$$

(34)

From the coefficient of $E_h$, which is $x_{2}a-x_1(a+{\mu }_h)=0$, this implies that; $x_2=x_1{\displaystyle \frac{(a+{\mu }_h)}{a}}$. Also, the coefficient of $I_m$, gives; $x_1{\rho }_h{\beta }_h{S}_h-x_4({\mu }_m+{\delta }_m)=0$, and implies that; $x_1=x_4{\displaystyle \frac{({\mu }_m+{\delta }_m)}{{\rho }_h{\beta }_h{S}_h}}$. The coefficient of $I_h$, have; $x_3{\rho }_m{\beta }_m{S}_m-x_2(b+{\mu }_h+{\delta }_h)=0$, putting $x_2$ in the coefficient of $I_h$, gives;

$x_3{\rho }_m{\beta }_m{S}_m-x_1{\displaystyle \frac{(a+{\delta }_m)(b+{\mu }_h+{\delta }_h)}{a}}=0$ and further simplification, resolved as;

$$x_3={\displaystyle \frac{\left[(a+{\delta }_m)(b+{\mu }_h+{\delta }_h)\right]x_1}{{\rho }_m{\beta }_m{S}_{m}a}}.$$

(35)

Therefore by substituting $x_1=x_4{\displaystyle \frac{({\mu }_m+{\delta }_m)}{{\rho }_h{\beta }_h{S}_h}}$ and equation (35) into equation (34), gives;

$${\displaystyle \frac{dV}{dt}}=\left[{\displaystyle \frac{x_4\theta -(a+{\delta }_m)(b+{\mu }_h+{\delta }_h)({\mu }_m+\theta )({\mu }_m+{\delta }_m)x_4}{{\rho }_m{\beta }_m{S}_m{\rho }_h{\beta }_h{S}_{h}a}}\right]E_m.$$

(36)

Further simplification gives;

$${\displaystyle \frac{dV}{dt}}=\left[{\displaystyle \frac{{\rho }_m{\beta }_m{S}_m{\rho }_h{\beta }_h{S}_{h}a\theta }{(a+{\delta }_m)(b+{\mu }_h+{\delta }_h)({\mu }_m+\theta )({\mu }_m+{\delta }_m)}}-1\right]x_4{E}_m.$$

(37)

Where;

$S_h={\displaystyle \frac{{\Lambda }_h}{{\mu }_h}}$ and $S_m={\displaystyle \frac{{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right)}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}.$

Suppose $x_4=1$, then;

$${\displaystyle \frac{dV}{dt}}=\left[R_0-1\right]E_m.$$

(38)

Thus, if $R_0\le 0$, then $\left[R_0-1\right]E_m\le 0$. Therefore, ${\displaystyle \frac{dV}{dt}}=0$. Furthermore, ${\displaystyle \frac{dV}{dt}}=0$ iff $E_m=0$, this shows that DFEP,

$E_{*}=\left({\displaystyle \frac{{\Lambda }_h}{{\mu }_h}},0,0,0,0,{\displaystyle \frac{{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right)}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}},0,0,{\displaystyle \frac{{\Lambda }_m{\rho }_q{\beta }_m}{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}\right)$ is globally asymptotically stable (GAS).

Remark 1.If $R_0<1$, the disease-free equilibrium is stable and the endemic equilibrium is absent, according to the stability analysis of the reproduction number.

Remark 2.In the event where $R_0>1$, the endemic equilibrium may exist and be stable, while the disease-free equilibrium is unstable. Hence the proved.

□

3. Sensitivity Analysis of the Model Using $R_0$

Sensitivity analysis of the model using the basic reproduction number, $R_0$, helps identify which parameters most significantly affect the spread of malaria. By analysing how variations in $R_0$ respond to changes in different parameters, researchers can pinpoint key areas for intervention. This information is crucial for optimizing control strategies and improving the effectiveness of insecticide use and community awareness programs.

Sensitivity analysis was used, following the methods described in [32], to ascertain the relative influence of each parameter on the spread of malaria. The normalised forward sensitivity index, as described in [33], is a variable $R_0$ that depends on the differentiable parameter $\gamma$. The analytical outcome of the sensitivity analysis of $R_0$ is obtained by computing ${\mathrm{X}}_{\gamma }^{R_0}={\displaystyle \frac{\partial R_0}{\partial \gamma }}\times {\displaystyle \frac{\gamma }{R_0}}$ to each of the parameters contributing to $R_0$. Thus, each of the fundamental parameters in $R_0$ has the following sensitivity;

$R_0=\sqrt{{\displaystyle \frac{{\Lambda }_m\left(\delta +{\mu }_m+{\delta }_a\right){\rho }_m{\beta }_{m}a{\rho }_h{\beta }_h{\Lambda }_h\theta }{{{\mu }_m}^2+\left(\delta +{\delta }_a\right){\mu }_m+\delta {\rho }_q{\beta }_m}}}.$

$$\begin{gathered} {\displaystyle \frac{\delta }{R_0}}\times {\displaystyle \frac{\partial }{\partial \delta }}{R}_0=-0.000007066908525;{\displaystyle \frac{{\mu }_m}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\mu }_m}}{R}_0=-0.4998702261; \\ {\displaystyle \frac{{\Lambda }_m}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\Lambda }_m}}{R}_0=+0.5;{\displaystyle \frac{{\delta }_a}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\delta }_a}}{R}_0=+0.00004368634364;{\displaystyle \frac{a}{R_0}}\times {\displaystyle \frac{\partial }{\partial a}}{R}_0=+0.5; \\ {\displaystyle \frac{\theta }{R_0}}\times {\displaystyle \frac{\partial }{\partial \theta }}{R}_0=+0.5;{\displaystyle \frac{{\Lambda }_h}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\Lambda }_h}}{R}_0=+0.5;{\displaystyle \frac{{\rho }_m}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\rho }_m}}{R}_0=+0.5; \\ {\displaystyle \frac{{\rho }_h}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\rho }_h}}{R}_0=+0.5;{\displaystyle \frac{{\beta }_h}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\beta }_h}}{R}_0=+0.5;{\displaystyle \frac{{\beta }_m}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\beta }_m}}{R}_0=+0.4998972084; \\ {\displaystyle \frac{{\rho }_q}{R_0}}\times {\displaystyle \frac{\partial }{\partial {\rho }_q}}{R}_0=-0.0001027913635 \end{gathered}$$.

3.1. Interpretation of the Sensitivity Index

The sensitivity index provides insight into which parameters have the most significant impact on the model outcomes. By analysing these indexes, researchers can identify which factors most strongly influence mosquito susceptibility to insecticides and the effectiveness of community awareness programs. This understanding helps in prioritising interventions and resources for more effective malaria control strategies. The sensitivity index in Table 3 for $R_0$ shows the expansion of mosquitoes in the community was significantly influenced by the parameters with positive index, which are ${\Lambda }_m,{\Lambda }_h,{\delta }_a,{\beta }_m,{\beta }_h,{\rho }_m,{\rho }_h,a$ and $\theta$. As their values increase, the burden on mosquitoes in the community is reduced by parameters with negative index $\mu ,\delta ,$ and ${\rho }_q$.

In order to combat disease in a community, the sensitivity technique showed that public health sectors and non-governmental organisations (NGOs) should use control tactics to decrease positives and improve control of negative index values. The following part examines an optimal control model in light of this observation in order to determine the best strategy for controlling diseases spread by mosquitoes.

4. Optimal Control Technique

The optimal control technique involves using mathematical models to determine the best strategies for managing mosquito populations and reducing malaria transmission. By adjusting variables such as the timing and quantity of insecticide application, this can minimize mosquito resistance while maximizing the effectiveness of control measures. This approach also considers community awareness programs to enhance public participation and improve overall outcomes in malaria prevention [32].

According to [20,34,35], the optimal control model is a powerful mathematical technique that may be used to complex situation analysis. According to [36], the control influences the dynamics of the system by entering the model equations for the ordinary differential equations. The goal is to modify control in order to maximise or minimise a particular functional aim. Similar to this, the problem requires a set of control variables, $u\left(t\right)\in U$, a collection of state variables, $y\left(t\right)\in Y$, and a cost function, $J\left[y\left(t\right),u\left(t\right)\right]$, for a time $t_0\le t\le t_f$, where $t_0$ and $t_f$ are the beginning and ending times, respectively.

The desired outcome is to determine a control $u\left(t\right)$ and the related state variables $x\left(t\right)$ in order to minimise a given objective function. The optimal control approaches in this study focuses as;

(i)

$u_1$: malathion, propoxur and permethrin chemical ingredients,

(ii)

$u_2$ : lasting nets insecticides (INI) and

(iii)

$u_3$ : traditional techniques.

The non-linear differential equations corresponds to the malaria transmission model diagram in Figure 1 and integrated optimal control technique, gives;

$$\left.\begin{array}{c}{\displaystyle \frac{d{S}_h}{dt}}={\Lambda }_h+\epsilon A_h-\left((1-u_1){\rho }_h{\beta }_h{I}_m+{\mu }_h\right)S_h,\hfill \\ \\ {\displaystyle \frac{d{A}_h}{dt}}=\alpha R_h-(\epsilon +{\mu }_h)A_h,\hfill \\ \\ {\displaystyle \frac{d{E}_h}{dt}}=(1-u_1){\rho }_h{\beta }_h{S}_h{I}_m-({\mu }_h+u_2+a)E_h,\hfill \\ \\ {\displaystyle \frac{d{I}_h}{dt}}=u_{2}a{E}_h-((1-u_3)+{\mu }_h+b+{\delta }_h)I_h,\hfill \\ \\ {\displaystyle \frac{d{R}_h}{dt}}=u_{3}b{I}_h-({\mu }_h+\alpha )R_h,\hfill \\ \\ {\displaystyle \frac{d{S}_m}{dt}}={\Lambda }_m+u_1\delta A_q-\left((1-u_1){\rho }_m{\beta }_m{I}_h+{\mu }_m\right)S_m,\hfill \\ \\ {\displaystyle \frac{d{E}_m}{dt}}=(1-u_1){\rho }_m{\beta }_m{S}_m{I}_h-({\mu }_m+u_1+\theta )E_m,\hfill \\ \\ {\displaystyle \frac{d{I}_m}{dt}}=u_1\theta E_m-({\mu }_m+{\delta }_m)I_m,\hfill \\ \\ {\displaystyle \frac{d{A}_q}{dt}}=(1-u_1){\rho }_q{\beta }_m(S_m+{\xi }_m{I}_m)-\left(u_1+\delta +{\mu }_m+{\delta }_a\right)A_q.\hfill \end{array}\right\}$$

(39)

with the initial conditions;

$S_h\left(0\right)>0,A_h\left(0\right)\ge 0,E_h\left(0\right)\ge 0,I_h\left(0\right)\ge 0,R_h\left(0\right)\ge 0,S_m\left(0\right)>0,E_m\left(0\right)\ge 0,I_m\left(0\right)\ge 0$, and $A\left(q\right)\left(0\right)\ge 0$.

In order to minimise the cost of control $u_1\left(t\right),u_2\left(t\right),and{u}_3\left(t\right)$ and to determine the optimal control values, equation (39) aims to minimise the total number of exposed individuals, infected individuals, and mosquitoes. $u^{*}=(u_1^{*},u_2^{*},u_3^{*}$ of the controls $u=(u_1,u_2,u_3)$ such that the associated state trajectories $(S_h^{•},A_h^{•},E_h^{•},I_h^{•},R_h^{•},S_m^{•},E_m^{•},I_m^{•},A_q^{•})$ are optimal control with initially conditions. The following is the objective function;

$$minJ(u_1,u_2,u_3)={\int }_{t_0}^{t_f}\left(C_1{E}_h\left(t\right)+C_2{I}_h\left(t\right)+C_3{Z}_m\left(t\right)+{\displaystyle \frac{1}{2}}\left(d_1{u}_1^2+d_2{u}_2^2+d_3{u}_3^2\right)\right)dt.$$

(40)

Subject to (41), where $t_f$ regards as final time and $Z_m=(E_m+I_m,)$ with $C_1,C_2,C_3,d_1,d_2$ and $d_3$ are positive weight constants and the choice of this study control agree with the idea in [37,38]. The $d_1{u}_1^2\left(t\right)+d_2{u}_2^2\left(t\right)+d_3{u}_3^2\left(t\right)$ are the costs associated with the use of malathion, propoxur and permethrin chemical ingredients, lasting nets insecticides (LNI), and traditional techniques. Consequently, the objective function makes it possible to maximise the Hamiltonian related to the optimum control problem. Afterwards, an optimal control $u_1^{*},u_2^{*}$, and $u_3^{*}$ satisfying equation (39) was obtained using the Pontryagin maximum principle [39];

$$J\left(u_1^{*},u_2^{*},u_3^{*}\right)=min\left\{J(u_1,u_2,u_3):u_1,u_2,u_3\in U\right\}.$$

(41)

Where $U=\left\{(u_1,u_2,u_3)\right\}$ such that $u_1,u_2,u_3$ are measurable with $0\le u_1\le 1,0\le u_2\le 1,0\le u_3\le 1,$ for $t\in \left[0,t_f\right]$ is the control set. The controls $u_1\left(t\right),u_2\left(t\right),u_3\left(t\right)$ In [40], are bdd Lebesgue integrable functions? $u_1\left(t\right),u_2\left(t\right),and{u}_3\left(t\right)$ are optimal controls that meet the requirements of Pontryagin’s maximal principle [41]. Equations (39) and (40) are transformed into a pointwise Hamiltonian H minimisation problem using this method, which yields the controls $u_1\left(t\right),u_2\left(t\right),and{u}_3\left(t\right)$. The Hamiltonian derived from the cost functional of equation (40) and the governing dynamics of equation (39) is used to determine the optimality conditions. Consequently, the following is the Hamiltonian H;

$$\left.\begin{array}{c}H=C_1{E}_h\left(t\right)+C_2{I}_h\left(t\right)+C_3{Z}_m\left(t\right)+{\displaystyle \frac{1}{2}}\left(d_1{u}_1^2+d_2{u}_2^2+d_3{u}_3^2\right)\hfill \\ +{\lambda }_1\left[{\Lambda }_h+\epsilon A_h-\left((1-u_1){\rho }_h{\beta }_h{I}_m+{\mu }_h\right)S_h\right]\hfill \\ +{\lambda }_2\left[\alpha R_h-(\epsilon +{\mu }_h)A_h\right]\hfill \\ +{\lambda }_3\left[(1-u_1){\rho }_h{\beta }_h{S}_h{I}_m-({\mu }_h+u_2+a)E_h\right]\hfill \\ +{\lambda }_4\left[u_{2}a{E}_h-((1-u_3)+{\mu }_h+b+{\delta }_h)I_h\right]\hfill \\ +{\lambda }_5\left[u_{3}b{I}_h-({\mu }_h+\alpha )R_h\right]\hfill \\ +{\lambda }_6\left[{\Lambda }_m+u_1\delta A_q-\left((1-u_1){\rho }_m{\beta }_m{I}_h+{\mu }_m\right)S_m\right]\hfill \\ +{\lambda }_7\left[(1-u_1){\rho }_m{\beta }_m{S}_m{I}_h-({\mu }_m+u_1+\theta )E_m\right]\hfill \\ +{\lambda }_8\left[u_1\theta E_m-({\mu }_m+{\delta }_m)I_m\right]\hfill \\ +{\lambda }_9\left[(1-u_1){\rho }_q{\beta }_m(S_m+{\xi }_m{I}_m)-\left(u_1+\delta +{\mu }_m+{\delta }_a\right)A_q\right]\hfill \end{array}\right\}.$$

(42)

Where ${\lambda }_i,i=1,2,2,...,9$ are the adjoints that correspond to the state variables. Taking into account, the equation (42) can be solve through appropriate partial derivatives of the Hamiltonian to the associated state variable.

Theorem 5.

Using the idea in [42,43]. Let the optimal control $u_1,u_2,u_3$ and $(S_h^{•},A_h^{•},E_h^{•},I_h^{•},R_h^{•},S_m^{•},E_m^{•},I_m^{•}$ and $A_q^{•})$ have a solution of the corresponding state variables that minimize J over U, and adjoint variables ${\lambda }_i,i=1,2,3,...9$ satisfying ${\displaystyle \frac{d{\lambda }_i}{dt}}=-{\displaystyle \frac{\partial H}{\partial i}}$.

Proof of Theorem 5.

According to [39], the integrand of the objective functional $u_1,u_2,u_3$ is a convex function of $u_1,u_2$ and $u_3$. The system satisfies the Lipshitz property with regard to the variables $u_1,u_2,u_3$ and $(S_h^{•},A_h^{•},E_h^{•},I_h^{•},R_h^{•},S_m^{•},E_m^{•},I_m^{•}$ and $A_q^{*})$ since the solution of (39) is bounded. Consequently, the optimal include $u_1^{*},u_2^{*},u_3^{*}$ exists. By differentiating the Hamiltonian function and evaluating it at the optimal control, the governing equations of the adjoint variables are derived. The adjoint system can therefore be expressed as;

$$\left.\begin{array}{c}{\displaystyle \frac{d{\lambda }_1}{dt}}=-{\displaystyle \frac{\partial H}{\partial S_h}}={\lambda }_1\left(\left(1-u_1\right){\rho }_h{\beta }_h{I}_m+{\mu }_h\right)-{\lambda }_3\left(1-u_1\right){\rho }_h{\beta }_h{I}_m,\hfill \\ \\ {\displaystyle \frac{d{\lambda }_2}{dt}}=-{\displaystyle \frac{\partial H}{\partial A_h}}={\lambda }_2\left(\epsilon +{\mu }_h\right)-{\lambda }_1ϵ,\hfill \\ \\ {\displaystyle \frac{d{\lambda }_3}{dt}}=-{\displaystyle \frac{\partial H}{\partial E_h}}=-C_1+{\lambda }_3\left({\mu }_h+u_2+a\right)-{\lambda }_4{u}_{2}a,\hfill \\ \\ {\displaystyle \frac{d{\lambda }_4}{dt}}=-{\displaystyle \frac{\partial H}{\partial I_h}}=-C_2+{\lambda }_4\left(1-u_3+{\mu }_h+b+{\delta }_h\right)-{\lambda }_5{u}_{3}b+{\lambda }_6\left(1-u_1\right){\rho }_m{\beta }_m{S}_m-{\lambda }_7\left(1-u_1\right){\rho }_m{\beta }_m{S}_m,\hfill \\ \\ {\displaystyle \frac{d{\lambda }_5}{dt}}=-{\displaystyle \frac{\partial H}{\partial R_h}}{\lambda }_5\left({\mu }_h+\alpha \right)-{\lambda }_2\alpha ,\hfill \\ \\ {\displaystyle \frac{d{\lambda }_6}{dt}}=-{\displaystyle \frac{\partial H}{\partial S_m}}={\lambda }_6\left(\left(1-u_1\right){\rho }_m{\beta }_m{I}_h+{\mu }_m\right)-{\lambda }_7\left(1-u_1\right){\rho }_m{\beta }_m{I}_h-{\lambda }_9\left(1-u_1\right){\rho }_q{\beta }_m,\hfill \\ \\ {\displaystyle \frac{{\lambda }_7}{dt}}=-{\displaystyle \frac{\partial H}{\partial E_m}}=-C_3+{\lambda }_7\left({\mu }_m+u_1+\theta \right)-{\lambda }_8{u}_1\theta ,\hfill \\ \\ {\displaystyle \frac{{\lambda }_8}{dt}}=-{\displaystyle \frac{\partial H}{\partial I_m}}=-C_3+{\lambda }_1\left(1-u_1\right){\rho }_h{\beta }_h{S}_h-{\lambda }_3\left(1-u_1\right){\rho }_h{\beta }_h{S}_h+{\lambda }_8\left({\mu }_m+{\delta }_m\right)-{\lambda }_9\left(1-u_1\right){\rho }_q{\beta }_m{\xi }_m,\hfill \\ \\ {\displaystyle \frac{{\lambda }_9}{dt}}=-{\displaystyle \frac{\partial H}{\partial A_q}}={\lambda }_9\left(u_1+\delta +{\mu }_m+{\delta }_a\right)-{\lambda }_6{u}_1\delta .\hfill \end{array}\right\}$$

(43)

The transversality conditions are given as;

$\lambda \left(t_f\right)$= 1,...,9.

Thus, the result give the control measure as set $(u_1,u_2,u_3)$, which follows with transversality conditions ${\lambda }_i\left(t_f\right)=1,2,...,9$, and the set of control as $(u_1,u_2,u_3)$, gives;

$\begin{array}{c}{u}_1=max\left\{0,min(1,{\zeta }_1),\right\}\hfill \\ u_2=max\left\{0,min(1,{\zeta }_2),\right\}\hfill \\ u_3=max\left\{0,min(1,{\zeta }_3).\right\}\hfill \end{array}$

Thus, the arguments, gives;

$$\begin{array}{c}{u}_1=\left\{\begin{array}{c}{\zeta }_1,if0<{\zeta }_1<1,\hfill \\ 0,if{\zeta }_1\le 0,\hfill \\ 1,if{\zeta }_1\ge 1,\hfill \end{array}\right\}\hfill \\ u_2=\left\{\begin{array}{c}{\zeta }_2,if0<{\zeta }_2<1,\hfill \\ 0,if{\zeta }_2\le 0,\hfill \\ 1,if{\zeta }_2\ge 1,\hfill \end{array}\right\}\hfill \\ u_3=\left\{\begin{array}{c}{\zeta }_3,if0<{\zeta }_3<1,\hfill \\ 0,if{\zeta }_3\le 0,\hfill \\ 1,if{\zeta }_3\ge 1.\hfill \end{array}\right\}\hfill \end{array}$$

(44)

Therefore;

$$\left.\begin{array}{c}{\zeta }_1={\displaystyle \frac{{\rho }_m{\beta }_m(I_m{\xi }_m+S_m){\lambda }_9-{\rho }_m{\beta }_m{I}_h{S}_m({\lambda }_6-{\lambda }_7)+{\lambda }_9{A}_q-{\rho }_h{\beta }_h{I}_m{S}_h({\lambda }_1-{\lambda }_3)-\delta {\lambda }_6{A}_q-E_m(\theta {\lambda }_8-{\lambda }_7)}{d_1}},\hfill \\ \\ {\zeta }_2={\displaystyle \frac{(a{\lambda }_4+{\lambda }_3)E_h}{d_2}},\hfill \\ \\ {\zeta }_3={\displaystyle \frac{(b{\lambda }_5+{\lambda }_4)I_h}{d_3}}.\hfill \end{array}\right\}$$

(45)

Furthermore, a sufficient condition for the Hamilton function’s second derivative with respect to $u_1,u_2,u_3$. This implies that ${\displaystyle \frac{{\partial }^{2}H}{\partial u_1^2}}=0\Rightarrow {\displaystyle \frac{{\partial }^{2}H}{\partial u_1^2}}=d_1,{\displaystyle \frac{{\partial }^{2}H}{\partial u_2^2}}=0\Rightarrow {\displaystyle \frac{{\partial }^{2}H}{\partial u_2^2}}=d_2$ and ${\displaystyle \frac{{\partial }^{2}H}{\partial u_3^2}}=0\Rightarrow {\displaystyle \frac{{\partial }^{2}H}{\partial u_3^2}}=d_3$ This indicates that the optimal problem is minimised at $u_1,u_2,u_3$.

Combining the established control set with the initial and transversality conditions allows one to predict the model equations with optimal control strategies in (39) and the adjoint variable as:

$$\left.\begin{array}{c}{\displaystyle \frac{d{S}_h}{dt}}={\Lambda }_h+\epsilon A_h-\left((1-u_1){\rho }_h{\beta }_h{I}_m+{\mu }_h\right)S_h,\hfill \\ {\displaystyle \frac{d{A}_h}{dt}}=\alpha R_h-(\epsilon +{\mu }_h)A_h,\hfill \\ {\displaystyle \frac{d{E}_h}{dt}}=(1-u_1){\rho }_h{\beta }_h{S}_h{I}_m-({\mu }_h+u_2+a)E_h,\hfill \\ {\displaystyle \frac{d{I}_h}{dt}}=u_{2}a{E}_h-((1-u_3)+{\mu }_h+b+{\delta }_h)I_h,\hfill \\ {\displaystyle \frac{d{R}_h}{dt}}=u_{3}b{I}_h-({\mu }_h+\alpha )R_h,\hfill \\ {\displaystyle \frac{d{S}_m}{dt}}={\Lambda }_m+u_1\delta A_q-\left((1-u_1){\rho }_m{\beta }_m{I}_h+{\mu }_m\right)S_m,\hfill \\ {\displaystyle \frac{d{E}_m}{dt}}=(1-u_1){\rho }_m{\beta }_m{S}_m{I}_h-({\mu }_m+u_1+\theta )E_m,\hfill \\ {\displaystyle \frac{d{I}_m}{dt}}=u_1\theta E_m-({\mu }_m+{\delta }_m)I_m,\hfill \\ {\displaystyle \frac{d{A}_q}{dt}}=(1-u_1){\rho }_q{\beta }_m(S_m+{\xi }_m{I}_m)-\left(u_1+\delta +{\mu }_m+{\delta }_a\right)A_q,\hfill \\ {\displaystyle \frac{d{\lambda }_1}{dt}}={\lambda }_1\left(\left(1-u_1\right){\rho }_h{\beta }_h{I}_m+{\mu }_h\right)-{\lambda }_3\left(1-u_1\right){\rho }_h{\beta }_h{I}_m,\hfill \\ {\displaystyle \frac{d{\lambda }_2}{dt}}={\lambda }_2\left(\epsilon +{\mu }_h\right)-{\lambda }_1ϵ,\hfill \\ {\displaystyle \frac{d{\lambda }_3}{dt}}=-C_1+{\lambda }_3\left({\mu }_h+u_2+a\right)-{\lambda }_4{u}_{2}a,\hfill \\ {\displaystyle \frac{d{\lambda }_4}{dt}}=-C_2+{\lambda }_4\left(1-u_3+{\mu }_h+b+{\delta }_h\right)-{\lambda }_5{u}_{3}b+{\lambda }_6\left(1-u_1\right){\rho }_m{\beta }_m{S}_m-{\lambda }_7\left(1-u_1\right){\rho }_m{\beta }_m{S}_m,\hfill \\ {\displaystyle \frac{d{\lambda }_5}{dt}}={\lambda }_5\left({\mu }_h+\alpha \right)-{\lambda }_2\alpha ,\hfill \\ {\displaystyle \frac{d{\lambda }_6}{dt}}={\lambda }_6\left(\left(1-u_1\right){\rho }_m{\beta }_m{I}_h+{\mu }_m\right)-{\lambda }_7\left(1-u_1\right){\rho }_m{\beta }_m{I}_h-{\lambda }_9\left(1-u_1\right){\rho }_q{\beta }_m,\hfill \\ {\displaystyle \frac{{\lambda }_7}{dt}}=-C_3+{\lambda }_7\left({\mu }_m+u_1+\theta \right)-{\lambda }_8{u}_1\theta ,\hfill \\ {\displaystyle \frac{{\lambda }_8}{dt}}=-C_3+{\lambda }_1\left(1-u_1\right){\rho }_h{\beta }_h{S}_h-{\lambda }_3\left(1-u_1\right){\rho }_h{\beta }_h{S}_h+{\lambda }_8\left({\mu }_m+{\delta }_m\right)-{\lambda }_9\left(1-u_1\right){\rho }_q{\beta }_m{\xi }_m,\hfill \\ {\displaystyle \frac{{\lambda }_9}{dt}}={\lambda }_9\left(u_1+\delta +{\mu }_m+{\delta }_a\right)-{\lambda }_6{u}_1\delta ,\hfill \\ {\zeta }_1={\displaystyle \frac{{\rho }_m{\beta }_m(I_m{\xi }_m+S_m){\lambda }_9-{\rho }_m{\beta }_m{I}_h{S}_m({\lambda }_6-{\lambda }_7)+{\lambda }_9{A}_q-{\rho }_h{\beta }_h{I}_m{S}_h({\lambda }_1-{\lambda }_3)-\delta {\lambda }_6{A}_q-E_m(\theta {\lambda }_8-{\lambda }_7)}{d_1}},\hfill \\ {\zeta }_2={\displaystyle \frac{(a{\lambda }_4+{\lambda }_3)E_h}{d_2}},\hfill \\ {\zeta }_3={\displaystyle \frac{(b{\lambda }_5+{\lambda }_4)I_h}{d_3}},\hfill \\ {\displaystyle \frac{\partial H}{\partial u_i}}=0,i=1,...3,\hfill \\ {\lambda }_1\left(t_f\right)={\lambda }_2\left(t_f\right)={\lambda }_3\left(t_f\right)={\lambda }_4\left(t_f\right)={\lambda }_5\left(t_f\right)={\lambda }_6\left(t_f\right)={\lambda }_7\left(t_f\right)={\lambda }_8\left(t_f\right)={\lambda }_9\left(t_f\right)=0.\hfill \end{array}\right\}$$

(46)

This complete the proved. □

5. Optimal Control Simulations of the Model

Optimal control simulations can be used to determine the most effective strategies for applying insecticides while minimising their use. By simulating different control scenarios, researchers can identify the best approach to reduce mosquito populations and interrupt malaria transmission. These simulations help in evaluating the costs, effectiveness, and potential resistance development. The system (1), optimal control technique in (39), adjoint system in (43), transversality boundary conditions, and initial conditions are $S_h\left(0\right)=20,000$,25], $E_h\left(0\right)=20$,40], $I_h\left(0\right)=50$,25], $A_h\left(0\right)=120,R_h\left(0\right)=100$,25], $A_q\left(0\right)=50,S_m\left(0\right)=1,000$,40], $E_m\left(0\right)=20,$ [40] and $I_m\left(0\right)=30$,40], while The cost weight constants corresponding to the state variables are assumed to be $C_1=C_2=C_3=1,d_1=40,d_2=0.5,d_3=45,$ and the set of adjoints ${\lambda }_i\left(t_f\right),i=1,...,9$. However, Table 4 contains the parameter values. Optimal control analysis is solved using numerical simulations using MATLABR2023a [40]. Using the fourth-order forward Runge-Kutta method, the equation (39) is solved given the initial values of state variables. The adjoint equations are then solved using the backward fourth-order Runge-Kutta method. [40,41,43]. According to [44] the following optimal strategies for control are used to assess how control measures affect the decreasing trends in malaria infections as;

(i)

Strategy A :malathion, propoxur and permethrin chemical ingredients $u_1$,

(ii)

Strategy B: lasting nets insecticide(LNI) $u_2$,

(iii)

Strategy C: traditional techniques $u_3$,

(iv)

Strategy D: combination of $u_1$ and $u_2$,

(v)

Strategy E: combination of $u_1$ and $u_3$,

(vi)

Strategy F: combination of $u_2$ and $u_3$,

(vii)

Strategy G: combination of $u_1,u_2$ and $u_3$.

6. Results and Discussion

The results indicate that implementing optimal insecticide application strategies can significantly reduce mosquito populations while minimising resistance development. Additionally, enhancing community awareness about malaria transmission dynamics leads to more effective control measures and increased participation in prevention efforts. The combination of these strategies creates a synergistic effect, improving overall public health outcomes in malaria-endemic regions. The sensitivity of mosquitoes strains to insecticides and traditional method, was examined. In Figure 2, the combined control result of the investigation showed that mosquitoes revealed no possible resistance to malathion, propoxur and permethrin as $u_1$ in this studied. It was also found that control $u_1$ showed no significant difference on susceptible unaware human, as shown in Figure 3. Subsequently, the result investigation revealed that $u_1$ and $u_3$ was found to be very effective toward mosquitoes in this study and agrees with the studies in [14] that found insecticide and awareness to be very active against the mosquito vector tested. The community awareness on mosquitoes is effective via campaign on mosquitoes to combat malaria, and can be deduced from Figure 4 and Figure 5 of the recovered individuals. Furthermore, the findings shows from Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15, Figure 16, Figure 17 and Figure 18 that employed a mathematical model to examine the resistance patterns of mosquitoes to commonly used insecticides, traditional method and evaluated the community’s awareness and practices related to malaria prevention proved to be effective. Based on the findings, it is recommended to implement targeted insecticide rotation strategies to manage resistance in mosquito populations effectively. Additionally, enhancing community education programs about malaria transmission through media campaign can increase awareness and encourage proactive participation in prevention efforts. Collaborative efforts between researchers, public health officials, and local communities are essential for sustainable control and reduction of malaria cases.

Strategy A: $u_1\ne 0,u_2=0,u_3=0$. The objective function J is optimised using controls $u_1\ne 0$, while $u_2$ and $u_3$ are set to zero in this method. In this case, the aim is to check the contribution of each combination intervention in reducing the mosquitoes and malaria transmission rate in a community. The results are shown in Figure 2, Figure 18 where the combination of $u_1$ (malathion, propoxur, and propoxur) is effective on susceptible,exposed, infected and aquatic mosquitoes, respectively. Also, from Figure 16, it can be deduced that only the insecticide $u_1$ is used on infected mosquitoes and does not work perfectly in reducing the number of vectors. This shows that relying on a single treatment only will not provide effective results. Subsequently, mass awareness and public campaigns on mosquito control are needed. Simply, this is because using insecticides to kill mosquitoes will help in reducing the susceptible. exposed, infected and aquatic mosquito population. As a result, the number of aquatic mosquitoes in the population drops, hence minimising the malaria transmission.

Strategy B: $u_2\ne 0,u_1-0,u_3=0$. The objective function J is optimized using this technique using the control $u_2\ne 0$, while $u_1=0,u_3=0$. The control combination strategies on exposed human in Figure 6 show no significant difference with that of Figure 7, where mosquitoes seem to be difficult to eliminate on the same dynamics in time. But in Figure 9, where $u_2\ne 0$ as shown, the combination control strategies revealed that lasting net insecticide effective as the number of exposed and infected mosquitoes are decreasing compared to that of aquatic mosquitoes in Figure 17 and Figure 18 due to the activities of humans to eradicate infected mosquitoes in a community through the use of insecticides. In this scenario, from Figure 14 and Figure 16, it can be deduced where mosquitoes seem to be eliminated using one or two interventions. Finally, it is recommended that these three intervention strategies be implemented simultaneously, especially in endemic areas, to effectively control malaria disease.

Strategy C: $u_3\ne 0,u_2=0,u_2=0$. Is optimised using the objective function J in this technique, whereas the other controls are set to zero in Figure 13. Also, from Figure 10, the adult mosquitoes seemed to have been eradicated, and the number of susceptible mosquitoes dramatically decreased throughout the first 100 days. Thus, when the single control strategy in Figure 8 with traditional method in Figure 9 is used to combat susceptible mosquitoes, the threat will be lessened compared to when control is not used.

Strategy D: $u_1\ne 0,u_2\ne 0,u_3=0$. Is optimised using the objective function J in this technique, while $u_3$ = 0 is set to zero in Figure 6. Also, from Figure 10, the number of recovered human seemed to have been progressed, and the number of exposed mosquitoes dramatically decreased throughout the first 100 days in Figure 14. Thus, when used lasting net insecticide $u_2$ and $u_3$ as combined control strategy in Figure 11, the threat will be lessened compared without control.

Strategy E: combination of $u_1\ne 0,u_3\ne 0$. The objective functional J is optimised using this technique, the control measures $(u_1,u_3)$, and $u_2=0$. The strategy’s applications’ outcomes are displayed in Figure 15 and Figure 17, it shows that the mosquito populations of exposed and aquatic are declining as compared to no control.

Strategy F: combination of $u_2\ne 0,u_3\ne 0$. The objective functional J is optimised using this technique, the control measures $(u_2,u_3)$, and $u_1=0$. The control applications are displayed in Figure 9 and Figure 11, it shows that the mosquito populations are declining as recovered individuals increasing.

Strategy G: combination of both A, B, and C with $u_1\ne 0,u_2\ne 0,u_3\ne 0$. The objective functional J is optimised using this technique, the control measures $(u_1,u_2)$, and $u_3$. The strategy’s applications outcomes are displayed in Figure 13 and Figure 17 it shows that the mosquito populations of susceptible, and aquatic are declining as compared to no control. Thus, the intervention technique is successful in reducing the mosquito population in time. Therefore, government decision-makers and all stakeholders are considered in applying all strategies to combat malaria in the specified time. Finally, there is a need to modify insecticide tactics and improve educational initiatives in order to combat the spread of malaria worldwide.

7. Conclusions

In conclusion, mathematical modelling and optimal control strategies play a crucial role in understanding and managing insecticide susceptibility in mosquitoes. By integrating community awareness campaigns about malaria transmission dynamics, this can enhance the effectiveness of control measures and reduce the incidence of malaria. Continued research and collaboration are essential to develop adaptive strategies that cater to evolving challenges in mosquito control and public health education. Finally, It is recommended to enhance community education and awareness programs about the proper use of insecticides to ensure effective mosquito control. Additionally, investing in research to develop new insecticides that mosquitoes have not yet developed resistance to can be beneficial. Collaboration with health organizations can also improve the implementation of optimal control techniques.