Introduction

JML

Journal of Mathematics Letters

2995-8075

Science Publications

10.31586/jml.2025.6104

JML-6104

Article

Mathematical modelling of the impact of HIV prevention strategies among female sex workers on public health in Burkina Faso

SOMDA

Serge M. A.

1 2 2 2 DABONÉ

Bernard E. A.

1 2 2 2 SANGARÉ

Boureima

1 2 2 2 * TRAORÉ

Sado

1 2 2 2

1 UFR SEA, Université Nazi-BONI, Burkina Faso

*Corresponding author at: UFR SEA, Université Nazi-BONI, Burkina Faso

26 06 2025

3 1 06 04 2025 23 05 2025 20 06 2025 26 06 2025

2025

This work is licensed under the Creative Commons Attribution International License (CC BY). http://creativecommons.org/licenses/by/4.0/

This article presents a mathematical model designed to simulate the impact of targeted interventions aimed at preventing HIV transmission among female sex workers (FSWs) and their clients, while also analyzing their effects on the health of the general population. The compartmental model distinguishes between high-risk populations (FSWs and their clients) and low-risk populations (sexually active men and women in the general population), and links prevention efforts in high-risk groups to the evolution of the epidemic in the low-risk population. The fundamental properties of the model, such as the positivity of solutions and the boundedness of the system, have been verified, and the basic reproduction number R<SUB>0</SUB> has been calculated. Finally, the stability of the model was studied using Varga’s theorem and the Lyapunov method. Simulation results show that targeted prevention among FSWs and their clients reduces HIV incidence in the general population. This framework provides a valuable tool for guiding policymakers in the design of effective strategies to combat the epidemic, especially relevant in the context of suspension of USAID funding.

Stability analysis; Key populations; HIV infection; Basic reproduction number; Mathematical modelling; Female sex workers

Introduction

Health systems are constantly evolving in order to improve the services offered to the population. When implementing these policies, it is important to be able to assess their impact on population’s health status. In the context of HIV, several prevention and treatment strategies have been defined to effectively control the spread of the disease [ 1]. According to UNAIDS, new HIV infections have been reduced by 59% since the peak in 1995 [ 2] and the current prevalence if 0.7% globally. However, the epidemic is most prevalent among key populations (KP), including female sex workers (FSW), men who have sex with men (MSM), injectable drug users (IDU), and prison mates [ 2]. As matter of fact, the prevalence of HIV infection was 3.6 times higher among sex workers, 10.7 times higher among MSMs, 7 times higher among IDUs, 14.7 times higher among transgender people and twice among the prison mates [ 3]. In Sub Saharan Africa, key populations accounted for 51% of new infections in 2021 [ 2]. Considering this new situation which represents a significant threat to disease control, a completely new response plan has been put in place by UNAIDS with KPs at the center of this new strategic plan [ 4].

Several mathematical models have been developed by researchers to provide a better understanding of the disease transmission in key population groups and to measure the effect of prevention and treatment strategies targeting them. In Kenya, for example, Omondi et al. [ 5] developed a mathematical model of HIV transmission between sex workers and injecting drug users to assess the effect of combination of Pre-Exposure Prophylaxis (PrEP) and antiretroviral therapy (ART) on the spread of HIV. In Cote d’Ivoire, Maheu-Giroux et al. [ 6] developed an age-stratified dynamic model of sexual and vertical transmission of HIV among the general population, FSWs and MSM. Their model was calibrated on detailed prevalence and intervention data (ART and condoms). Geidelberg et al. [ 7] looked at the role of PrEP on the spread of HIV among FSWs in Cotonou (Benin Republic). The authors also used a compartmental model involving PrEP and ART in high-risk (FSWs and clients) and low-risk populations. In Burkina Faso, Low et al. [ 8] formulated an initial mathematical model to assess the impact of combining ARTs and condoms on FSWs’ health status.

So far, there is no standardised approach to directly visualising how interventions targeting key populations can affect the general population that does not interact directly with them. The aim of the model proposed in this study is therefore to describe and simulate the potential benefits of HIV control programs for key populations, taking into account the indirect effects on the general population. This research extends the initial study by Low et al. [ 8]. To illustrate our approach at the national level, we examine the specific situation of Burkina Faso, a landlocked country in West Africa where HIV prevalence varies widely among demographic groups and identified key populations. In 2019, the national census reported a population of 20.5 million, currently estimated to be close to 23.5 million [ 9]. HIV prevalence in the general population is estimated at 0.6% [ 10], with specific rates of 0.8% for women and 0.5% for men aged 15-49 years. Among key populations, prevalence is estimated at 5.4% among sex workers, 1.9% among men who have sex with men, 1.3% among prisoners and 1.0% among injecting drug users [ 11].

The work in this paper is organized as followed: in Section 2, the model describing the transmission of HIV between high and low risk individuals is formulated. In Section 3 the mathematical analysis of the model is given. Application with numerical simulation and the model simulation results are presented in Sections 4. Finally, the discussion and conclusion are presented in Section 5.

Model formulation

In this section, a mathematical model is formulated. It describes the HIV transmission between two different risk groups, namely sex workers and their clients considered to be at high risk of infection and the general sexually active population considered to be at low risk of infection. This model is an extension of a model studied in [ 12] and is based on the modelling approach given in [ 5,13]. The transition between these risk groups is assumed to exist in this model. The model has twelve compartments. S, I, and A represent the level of infection according to a classical SI model. These are respectively the susceptible, the infectious and the AIDS class in the model. The type of group ${s w, c, m, f}$ represents the category of the participants which are respectively the FSWs, their clients, the general male population and the general female population. The infection rate $λ$ is determined from the average annual number of sexual partners ( $ψ$ ), the probability of HIV transmission during vaginal sexual intercourse ( $p_{M}$ for male and $p_{F}$ for female) with a partner in infectious (I) or AIDS (A) state. The parameter $Γ$ will permit to modify the risks of HIV sexual transmissions. Infection rates in at-risk populations are given by

λ_{s w} = p_{M} (1 - π) \frac{Γ_{1} ψ_{s w}^{I_{c}} I_{c} + Γ_{2} ψ_{s w}^{A_{c}} A_{c}}{N_{c}},

λ_{s w} = p_{F} [(1 - π) \frac{Γ_{1} ψ_{c}^{I_{s w}} I_{s w} + Γ_{2} ψ_{c}^{A_{s w}} A_{s w}}{N_{s w}} + \frac{Γ_{1} ψ_{c}^{I_{f}} I_{f} + Γ_{2} ψ_{c}^{A_{f}} A_{f}}{N_{f}}],

λ_{f} = p_{M} [\frac{Γ_{1} ψ_{f}^{I_{c}} I_{c} + Γ_{2} ψ_{f}^{A_{c}} A_{c}}{N_{c}} + \frac{Γ_{1} ψ_{f}^{I_{m}} I_{m} + Γ_{2} ψ_{f}^{A_{m}} A_{m}}{N_{m}}],

λ_{s w} = p_{F} \frac{Γ_{1} ψ_{m}^{I_{f}} I_{f} + Γ_{2} ψ_{m}^{A_{f}} A_{f}}{N_{f}} .

The parameter $π$ is added to consider the effect of preventive intervention targeting FSWs and their clients. The model assumes that all risk groups grow at the same rate over each period. The population sizes at time t for the risk groups are:

N_{s w} (t) = S_{s w} (t) + I_{s w} (t) + A_{s w} (t); N_{c} (t) = S_{c} (t) + I_{c} (t) + A_{c} (t);

N_{m} (t) = S_{m} (t) + I_{m} (t) + A_{m} (t); N_{f} (t) = S_{f} (t) + I_{f} (t) + A_{f} (t);

It is important to note that the intervention is implemented in the group of high-risk individuals. The full list of model parameters is summarized inTable table 1.

The flow and interactions between the compartments are described inFigure figure 1. In this figure, black arrows indicate the evolution of individuals in the course of the infection while red arrows indicate sexual infection across groups.

Table 1

Table 1. Description of parameters used in the model.

Symbol	Description
d	AIDS fatality rates

Figure 1

A compartmental representation of the model for HIV transmission.

The above description gives the following non-linear ordinary differential equations: For sex workers and their clients

\{\begin{matrix} \frac{d S_{s w}}{d t} = Λ_{s w} - (λ_{s w} - μ) S_{s w}, \\ \frac{d I_{s w}}{d t} = λ_{s w} S_{s w} - (γ + μ) I_{s w} \\ \frac{d A_{s w}}{d t} = γ I_{s w} - (μ + d) A_{s w}, \\ \frac{d S_{c}}{d t} = Λ_{s w} - (λ_{s w} - μ) S_{s w}, \\ \frac{d I_{c}}{d t} = λ_{c} S_{c} - (γ + μ) I_{c}, \\ \frac{d A_{c}}{d t} = γ I_{c} - (μ + d) A_{c} . \end{matrix}

(1)

For men and women in general population

\{\begin{matrix} \frac{d f}{d t} = Λ_{f} - (λ_{f} - μ) S_{f}, \\ \frac{d I_{f}}{d t} = λ_{f} S_{f} - (γ + μ) I_{f}, \\ \frac{d A_{f}}{d t} = γ I_{f} - (μ + d) A_{f}, \\ \frac{d S_{m}}{d t} = Λ_{m} - (λ_{m} - μ) S_{m}, \\ \frac{d I_{m}}{d t} = λ_{m} S_{m} - (γ + μ) I_{m}, \\ \frac{d A_{m}}{d t} = γ I_{m} - (μ + d) A_{m} . \end{matrix}

(2)

The system of equations in (1) — (2) is subject to the following initial conditions $, S_{s w} (0) \geq 0, I_{s w} ((0) \geq 0, A_{s w} ((0) \geq 0, S_{c} (0) \geq 0, I_{c} (0) \geq 0, A_{c} (0) \geq 0, S_{m} (0) \geq 0, I_{m} (0)) \geq 0, A_{m} (0) \geq 0, S_{f} (0) \geq 0, I f (0) \geq 0$ and $A_{f} (0)$ ≥ 0.

Model analysis3.1. Positivity and boundedness properties

In this subsection, we determine the equilibrium points, and the basic reproduction number associated with System (1) — (2).

Lemma 3.1. If $S_{s w} (0)$ , $I_{s w} (0)$ , $A_{s w} (0)$ , $S_{c} (0)$ , $I_{c} (0)$ , $A_{c} (0)$ , $S_{m} (0)$ , $I_{m} (0)$ , $A_{m} (0)$ , $S_{f} (0)$ , $I_{f} (0)$ , and $A_{f} (0)$ are non negative, then so are $S_{s w} (t),$ $I_{s w} (t)$ , $A_{s w} (t)$ , $S_{c} (t)$ , $I_{c} (t)$ , $A_{c} (t)$ , $S_{m} (t)$ , $I_{m} (t)$ , $A_{m} (t)$ , $S_{f} (t)$ , $I_{f} (t)$ , and $A_{f} (t)$ for all time $t > 0$ . Moreover, ${l i m}_{t \to + \infty} S u p N (t) \leq \frac{Λ}{μ}$ .

Proof 3.1. For that, we first use the contradiction that the state variable $S$ is nonnegative for all $t \leq 0$ . Let $e (t) = m i n {S_{s w} (t), I_{s w} (t), A_{s w} (t), S_{c O} (t), I_{c} (t), A_{c} (t), S_{m} (t), I_{m} (t), A_{m} (t), S_{f} (t), I_{f} (t), A_{f} (t)}$ and let us suppose that there is at least one $t_{1} > 0$ such that $e (t_{1}) = 0$ and $e (t_{1}) > 0$ for all $t \in (0, t_{1})$ . Therefore, if $e (t) = S_{s w} (t)$ then each state of model (1) $&$ (2), is positive and from the first equation of (1) $&$ (2), we have

{\dot{S}}_{s w} (t) > μ S_{s w} (t)

(3)

It then follows that

0 = S_{s w} (t_{1}) > S_{s w} (0) e^{- μ t_{1}}

(4)

which leads to a contradiction.

Similar proof can be given for the other state variables. Thus, any solution of system (1) — (2) is nonnegative for $t \geq 0$ . Moreover, the total number of the population N(t) at any time is governed by

\dot{N} (t) = Λ — µ N (t) — d (A s w (t) + A c (t) + A f (t) + A m (t)) \leq Λ — µ N (t) .

(5)

Where $Λ = Λ_{s w} + Λ_{c} + Λ_{f} + Λ_{m}$ .

Thus, for the initial conditions $0 \leq N (0) \leq \frac{Λ}{μ}$ , by using Gronwall inequality [ 14,15], we get

0 \leq N (0) \leq \frac{Λ}{μ}

(6)

Hence, system (1) — (2) defines a dynamical on

Δ = \{(S_{s w} (t), I_{s w} (t), A_{s w} (t), S_{c} (t), I_{c} (t), A_{c} (t), S_{m} (t), I_{m} (t), A_{m} (t), S_{f} (t), I_{f} (t), A_{f} (t)) \in R_{+}^{12} : 0 \leq N (0) \leq \frac{Λ}{μ}\} .

3.2. Equilibrium points and basic reproduction number

In this subsection, we determine the equilibrium points, and the basic reproduction number associated with System (1) — (2).

3.2.1. Equilibrium points

Let ε = (S_sw, I_sw, A_sw, S_c, I_c, A_c, S_f , I_f , A_f , S_m, I_m, A_m) be the equilibrium point of System. This corresponds to the values where the equations (1) — (2) are all null.

Let ε₀ and ε^∗ be respectively the disease-free and endemic equilibrium points of System (1) — (2). The disease-free equilibrium corresponds to the case were there is not infected individual. In this case, we have $I_{s w} = 0, A_{s w} = 0, I_{c} = 0, A_{c} = 0, I_{f} = 0, A_{f} = 0, I_{m} = 0, A_{m} = 0$ and the disease-free equilibrium given by:

ε_{0} = (\frac{Λ_{s w}}{μ}, 0,0, \frac{Λ_{c}}{μ}, 0,0, \frac{Λ_{f}}{μ}, 0,0, \frac{Λ_{m}}{μ}, 0,0)

(7)

We design by $I_{s w}^{*}, A_{s w}^{*}, I_{c}^{*}, A_{c}^{*}, I_{f}^{*}, A_{f}^{*}, I_{m}^{*}$ and $A_{s w}^{*}$ the infectious at endemic equilibrium points, so

ε * = (S_{s w}^{*}, I_{s w}^{*}, A_{s w}^{*}, S_{c}^{*}, I_{c}^{*}, A_{c}^{*}, S_{f}^{*}, I_{f}^{*}, A_{f}^{*}, S_{m}^{*}, I_{m}^{*}, A_{m}^{*})

(8)

With

S_{s w}^{*} = \frac{Λ_{s w}}{μ} - \frac{(μ + γ) (μ + d)}{γ μ} A_{s w}^{*},

I_{s w}^{*} = \frac{μ + d}{γ} A_{s w}^{*},

S_{c}^{*} = \frac{Λ_{c}}{μ} - \frac{(μ + γ) (μ + d)}{γ μ} A_{c}^{*},

I_{c}^{*} = \frac{μ + d}{γ} A_{c}^{*},

S_{f}^{*} = \frac{Λ_{f}}{μ} - \frac{(μ + γ) (μ + d)}{γ μ} A_{f}^{*},

I_{f}^{*} = \frac{μ + d}{γ} A_{f}^{*},

S_{m}^{*} = \frac{Λ_{m}}{μ} - \frac{(μ + γ) (μ + d)}{γ μ} A_{m}^{*},

I_{m}^{*} = \frac{μ + d}{γ} A_{m}^{*} .

By putting the right hand side of system (1) — (2) to zero, and keeping each state variable different from zero ( $S_{s w} \neq 0$ , $I_{s w} \neq 0$ , $S_{c} \neq 0$ , $I_{c} \neq 0$ , $I_{f} \neq 0$ , $I_{f} \neq 0$ , $S_{m} \neq 0$ and $I_{m} \neq 0$ ).

We also make this assumption:

(H) : \forall I_{s w}^{*}, I_{c}^{*}, I_{f}^{*}, I_{m}^{*} \in R_{+}, (\begin{matrix} I_{s w} \\ I_{c} \\ I_{f} \\ I_{m} \end{matrix}) \leq (\begin{matrix} I_{s w}^{*} \\ I_{c}^{*} \\ I_{f}^{*} \\ I_{m}^{*} \end{matrix}) .

3.2.2. The Basic reproduction numbers

The threshold parameter R₀ is defined as the average cases of secondary infections generated by a single infectious individual in a completely susceptible population during his/her period of infectiousness [ 16,17]. In our model the infected classes correspond to states I_sw, A_sw, I_c, A_c, I_f, A_f, I_m and A_m. Thus, we can rewrite system (1) — (2) as

\dot{χ} = F (χ) - V (χ),

(9)

where $χ = (I_{s w}, A_{s w}, I_{c}, A_{c}, I_{f}, A_{f}, I_{m}, A_{m})$ ; $F$ is the rate of appearance of new infections in each class, and $V$ is the rate of transfer of individuals out of (for positive values) or into (for negative values) compartment by all other means. Hence, in our case, we have that

F = (\begin{matrix} [p_{M} (1 - π) \times \frac{Γ_{1} ψ_{s w}^{I_{c}} I_{c} + Γ_{2} ψ_{s w}^{A_{c}} A_{c}}{N_{c}}] \times S_{s w} \\ 0 \\ p_{F} [(1 - π) \times \frac{Γ_{1} ψ_{c}^{I_{s w}} I_{s w} + Γ_{2} ψ_{c}^{A_{s w}} A_{s w}}{N_{s w}} + \frac{Γ_{1} ψ_{c}^{I_{f}} I_{f} + Γ_{2} ψ_{c}^{A_{f}} A_{f}}{N_{f}}] \times S_{c} \\ 0 \\ p_{M} [\frac{Γ_{1} ψ_{f}^{I_{c}} I_{c} + Γ_{2} ψ_{f}^{A_{c}} A_{c}}{N_{c}} + \frac{Γ_{1} ψ_{f}^{I_{m}} I_{m} + Γ_{2} ψ_{f}^{A_{m}} A_{m}}{N_{m}}] \times S_{f} \\ 0 \\ p_{F} \times [\frac{Γ_{1} ψ_{m}^{I_{f}} I_{f} + Γ_{2} ψ_{m}^{A_{f}} A_{f}}{N_{f}}] \times S_{m} \\ 0 \end{matrix}),

and

V = (\begin{matrix} (μ + γ) I_{s w} \\ - γ I_{s w} + (μ + d) A_{s w} \\ (μ + γ) I_{c} \\ - γ I_{c} + (μ + d) A_{c} \\ (μ + γ) I_{f} \\ - γ I_{f} + (μ + d) A_{f} \\ (μ + γ) I_{m} \\ - γ I_{m} + (μ + d) A_{m} \end{matrix}) .

The jacobian matrix of $F$ and $V$ at the disease-free equilibrium $ε_{0}$ are given by

F = (\begin{matrix} 0 & 0 & F_{1,3} & F_{1,4} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ F_{3,1} & F_{3,2} & 0 & 0 & F_{3,5} & F_{3,6} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & F_{5,3} & F_{5,4} & 0 & 0 & F_{5,7} & F_{5,8} \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & F_{7,5} & F_{7,6} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{matrix}),

where:

F_{1,3} = p_{M} \times (1 - π) \times (Γ_{1} ψ_{s w}^{I_{c}}) \times \frac{S_{s w}^{0}}{S_{c}^{0}}; F_{1,4} = p_{M} \times (1 - π) \times (Γ_{2} ψ_{s w}^{A_{c}}) \times \frac{S_{s w}^{0}}{S_{c}^{0}};

F_{3,1} = p_{F} \times (1 - π) \times (Γ_{1} ψ_{c}^{I_{s w}}) \times \frac{S_{c}^{0}}{S_{s w}^{0}}; F_{3,2} = p_{F} \times (1 - π) \times (Γ_{2} ψ_{c}^{A_{s w}}) \times \frac{S_{c}^{0}}{S_{s w}^{0}};

F_{3,5} = p_{F} \times (Γ_{1} ψ_{c}^{I_{f}}) \times \frac{S_{c}^{0}}{S_{f}^{0}}; F_{3,6} = p_{F} \times (Γ_{2} ψ_{c}^{A_{f}}) \times \frac{S_{c}^{0}}{S_{f}^{0}};

F_{5,3} = p_{M} \times (1 - π) \times (Γ_{1} ψ_{f}^{I_{c}}) \times \frac{S_{f}^{0}}{S_{c}^{0}}; F_{5,4} = p_{M} \times (1 - π) \times (Γ_{2} ψ_{f}^{A_{c}}) \times \frac{S_{f}^{0}}{S_{c}^{0}};

F_{5,7} = p_{M} \times (1 - π) \times (Γ_{1} ψ_{f}^{I_{m}}) \times \frac{S_{f}^{0}}{S_{m}^{0}}; F_{5,8} = p_{M} \times (1 - π) \times (Γ_{2} ψ_{f}^{A_{m}}) \times \frac{S_{f}^{0}}{S_{m}^{0}};

F_{7,5} = p_{F} \times (Γ_{1} ψ_{m}^{I_{f}}) \times \frac{S_{m}^{0}}{S_{f}^{0}}; F_{7,6} = p_{F} \times (Γ_{2} ψ_{m}^{A_{f}}) \times \frac{S_{m}^{0}}{S_{f}^{0}} .

and

V = (\begin{matrix} (μ + γ) & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ - γ & (μ + d) & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & (μ + γ) & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & - γ & (μ + d) & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & (μ + γ) & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & - γ & (μ + d) & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & (μ + γ) & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & - γ & (μ + d) \end{matrix}),

respectively.

The matrix $V$ is invertible (non-zero determinant) and its inverse $V^{- 1}$ is defined by:

V^{- 1} = (\begin{matrix} \frac{1}{(μ + γ)} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{γ}{(μ + γ) (μ + d)} & \frac{1}{(μ + d)} & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{1}{(μ + γ)} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{γ}{(μ + γ) (μ + d)} & \frac{1}{(μ + d)} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{(μ + γ)} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{γ}{(μ + γ) (μ + d)} & \frac{1}{(μ + d)} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & \frac{1}{(μ + γ)} & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & \frac{γ}{(μ + γ) (μ + d)} & \frac{1}{(μ + d)} \end{matrix}),

So, the next-generation matrix $F V^{- 1}$ is given by:

(\begin{matrix} 0 & 0 & \frac{γ F_{1,4}}{(μ + γ) (μ + d)} + \frac{F_{1,3}}{(μ + γ)} & \frac{F_{1,4}}{(μ + d)} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{γ F_{3,2}}{(μ + γ) (μ +)} + \frac{F_{3,1}}{(μ + γ)} & \frac{F_{3,2}}{(μ + d)} & 0 & 0 & \frac{γ F_{3,6}}{(μ + γ) (μ + d)} + \frac{F_{3,5}}{(μ + γ)} & \frac{F_{3,6}}{(μ + d)} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{γ F_{5,4}}{(μ + γ) (μ + d)} + \frac{F_{5,3}}{(μ + γ)} & \frac{F_{5,4}}{(μ + d)} & 0 & 0 & \frac{γ F_{5,8}}{(μ + γ) (μ + d)} + \frac{F_{5,7}}{(μ + γ)} & \frac{F_{5,8}}{(μ + d)} \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{γ F_{7,6}}{(μ + γ) (μ + d)} + \frac{F_{7,5}}{(μ + γ)} & \frac{F_{7,6}}{(μ + d)} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{matrix}),

Let’s put

A = \frac{γ F_{1,4}}{(μ + γ) (μ + d)} + \frac{F_{1,3}}{(μ + γ)}; B = \frac{F_{1,4}}{(μ + d)}; C = \frac{γ F_{3,2}}{(μ + γ) (μ + d)} + \frac{F_{3,1}}{(μ + γ)}; D = \frac{F_{3,2}}{(μ + d)};

E = \frac{γ F_{3,6}}{(μ + γ) (μ + d)} + \frac{F_{3,5}}{(μ + γ)}; F = \frac{F_{3,6}}{(μ + d)}; G = \frac{γ F_{5,4}}{(μ + γ) (μ + d)} + \frac{F_{5,3}}{(μ + γ)}; H = \frac{F_{5,4}}{(μ + d)};

I = \frac{γ F_{5,8}}{(μ + γ) (μ + d)}; J = \frac{F_{5,8}}{(μ + d)}; K = \frac{γ F_{7,6}}{(μ + γ) (μ + d)} + \frac{F_{7,5}}{(μ + γ)}; L = \frac{F_{7,6}}{(μ + d)} .

The basic reproduction ratio is given by $R_{0} = ρ (- F V^{- 1})$ , the spectral radius of the next-generation matrix $F V^{- 1}$ . In our case, we have that

R_{0} = \{\begin{matrix} \frac{- β + \sqrt{β^{2} - 4 α}}{2}, i f β^{2} - 4 α > 0 \\ \frac{β}{2}, i f β^{2} - 4 α = 0 \end{matrix}

Where:

α = A C I K a n d β = - (I K + E G + A C)

3.3. Stability of equilibrium points

In this subsection, we prove the global stability of disease-free equilibrium $ε_{0}$ , when $R_{0} < 1$ and the global stability of endemic equilibrium $ε^{*}$ , when $R_{0} > 1$ .

Lemma 3.2. The disease-free equilibrium $ε_{0}$ of system (1) — (2) is globally asymptotically stable, when $R_{0} < 1$ and unstable if $R_{0} > 1$ .

Proof 3.2. Let us consider the infected classes $χ = (I_{s w}, A_{s w}, I_{c}, A_{c}, I_{f}, A_{f}, I_{m}, A_{m})$ . By the equations corresponding to these states, we have the linearization system at $ε_{0}$ given by:

\{\begin{array}{l} \frac{d I_{s w}}{d t} = p_{M} (1 - π) (\frac{Γ_{1} ψ_{s w}^{I_{c}} Λ_{s w}}{μ N_{c}}) I_{c} + p_{M} (1 - π) (\frac{Γ_{2} ψ_{s w}^{A_{c}} Λ_{s w}}{μ N_{c}}) A_{c} - (μ + γ) I_{s w}, \\ \frac{d A_{s w}}{d t} = γ I_{s w} - (μ + d) A_{s w}, \\ \frac{d I_{c}}{d t} = p_{F} (1 - π) (\frac{Γ_{1} ψ_{c}^{I_{s w}} Λ_{c}}{μ N_{s w}}) I_{s w} + p_{M} (1 - π) (\frac{Γ_{2} ψ_{c}^{A_{s w}} Λ_{c}}{μ N_{s w}}) A_{s w} + p_{F} (\frac{Γ_{1} ψ_{c}^{I_{f}} Λ_{c}}{μ N_{f}}) I_{f} + p_{F} (\frac{Γ_{2} ψ_{c}^{A_{f}} Λ_{c}}{μ N_{f}}) A_{f} - (μ + γ) I_{c}, \\ \frac{d A_{c}}{d t} = γ I_{c} - (μ + d) A_{c}, \\ \frac{d I_{f}}{d t} = p_{M} (\frac{Γ_{1} ψ_{f}^{I_{c}} Λ_{f}}{μ N_{c}}) I_{c} + p_{F} (\frac{Γ_{2} ψ_{f}^{A_{c}} Λ_{f}}{μ N_{c}}) A_{c} + p_{F} (\frac{Γ_{1} ψ_{f}^{I_{m}} Λ_{f}}{μ N_{m}}) I_{m} + p_{F} (\frac{Γ_{2} ψ_{f}^{A_{m}} Λ_{f}}{μ N_{m}}) A_{m} - (μ + γ) I_{f}, \\ \frac{d A_{f}}{d t} = γ I_{f} - (μ + d) A_{f}, \\ \frac{d I_{m}}{d t} = p_{F} (\frac{Γ_{1} ψ_{m}^{I_{f}} Λ_{m}}{μ N_{f}}) I_{f} + p_{F} (\frac{Γ_{2} ψ_{m}^{A_{f}} Λ_{m}}{μ N_{f}}) A_{f} - (μ + γ) I_{m}, \\ \frac{d A_{m}}{d t} = γ I_{m} - (μ + d) A_{m} . \end{array}

(10)

the matrix associated to the linearized system (11) is given by:

A = (\begin{matrix} - (μ + γ) & 0 & A_{1,3} & A_{1,4} & 0 & 0 & 0 & 0 \\ γ & - (μ + d) & 0 & 0 & 0 & 0 & 0 & 0 \\ A_{3,1} & A_{3,2} & - (μ + γ) & 0 & A_{3,5} & A_{3,6} & 0 & 0 \\ 0 & 0 & γ & - (μ + d) & 0 & 0 & 0 & 0 \\ 0 & 0 & A_{5,3} & A_{5,4} & - (μ + γ) & 0 & A_{5,7} & A_{5,8} \\ 0 & 0 & 0 & 0 & γ & - (μ + d) & 0 & 0 \\ 0 & 0 & 0 & 0 & A_{7,5} & A_{7,6} & - (μ + γ) & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & γ & - (μ + d) \end{matrix})

With

A_{1,3} = p_{M} (1 - π) (\frac{Γ_{1} ψ_{s w}^{I_{c}} Λ_{s w}}{μ N_{c}}) A_{1,4} = p_{M} (1 - π) (\frac{Γ_{2} ψ_{s w}^{A_{c}} Λ_{s w}}{μ N_{c}}) .

A_{3,1} = (1 - π) (\frac{Γ_{1} ψ_{c}^{I_{s w}} Λ_{c}}{μ N_{s w}}); A_{3,2} = p_{M} (1 - π) (\frac{Γ_{2} ψ_{c}^{A_{s w}} Λ_{c}}{μ N_{s w}}); A_{3,5} = p_{F} (\frac{Γ_{1} ψ_{c}^{I_{f}} Λ_{c}}{μ N_{f}}); A_{3,6} = p_{F} (\frac{Γ_{2} ψ_{c}^{A_{f}} Λ_{c}}{μ N_{f}}) .

A_{5,3} = p_{M} (1 - π) (\frac{Γ_{1} ψ_{f}^{I_{c}} Λ_{f}}{μ N_{c}}); A_{5,4} = p_{F} (1 - π) (\frac{Γ_{2} ψ_{f}^{A_{c}} Λ_{f}}{μ N_{c}}); A_{5,7} = p_{F} (\frac{Γ_{1} ψ_{f}^{I_{m}} Λ_{f}}{μ N_{m}}); A_{5,8} = p_{F} (\frac{Γ_{2} ψ_{f}^{A_{m}} Λ_{f}}{μ N_{m}}) .

A_{7,5} = p_{F} (\frac{Γ_{1} ψ_{m}^{I_{f}} Λ_{m}}{μ N_{f}}); A_{7,5} = p_{F} (\frac{Γ_{2} ψ_{m}^{A_{f}} Λ_{m}}{μ N_{f}}) .

It is important to note that $A$ has all off diagonal entries non-negative. This implies that $A$ is a Metzler matrix. The matrix A can be decomposed as follows $A = F_{1} + V_{1}$ , with:

F_{1} = (\begin{matrix} 0 & 0 & A_{1,3} & A_{1,4} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ A_{3,1} & A_{3,2} & 0 & 0 & A_{3,5} & A_{3,6} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & A_{5,3} & A_{5,4} & 0 & 0 & A_{5,7} & A_{5,8} \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & A_{7,5} & A_{7,6} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & - 0 \end{matrix});

and

V_{1} = (\begin{matrix} - (μ + γ) & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ γ & - (μ + d) & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & - (μ + γ) & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & γ & - (μ + d) & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & - (μ + γ) & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & γ & - (μ + d) & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & - (μ + γ) & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & γ & - (μ + d) \end{matrix}) .

The matrix $V_{1}$ is an invertible matrix and it invert $V_{1}^{- 1}$ is given by the matrix $V^{- 1}$ . We can also see that $F_{1} \geq 0$ and $V_{1}^{- 1} \geq 0$ because of $V_{1}^{- 1}$ is also a Metzer matrix. Thus, $R_{0} = ρ (- F_{1} V_{1}^{- 1})$ and from the Varga theorem [ 18]. the matrix $A$ is asymptotically stable. Therefore, disease free equilibrium $ε_{0}$ is globally asymptotically stable when $R_{0} < 0$ .

Lemma 3.3. The endemic equilibrium ε ∗ of system (1) − (2) is globally asymptotically stable, when $R_{0} > 0$ .

Proof 3.3. For this proof, we will use the same approach developed in [13, 15, 19-22].

Let $ε^{*} = (S_{s w}^{*}, I_{s w}^{*}, A_{s w}^{*}, S_{c}^{*}, I_{c}^{*}, A_{c}^{*}, S_{f}^{*}, I_{f}^{*}, A_{f}^{*}, S_{m}^{*}, I_{m}^{*}, A_{m}^{*})$ be the endemic equilibrium of systems (1) - (2).From systems (1) -(2), we have:

\begin{array}{r} \{\begin{array}{l} Λ_{s w} & = λ_{s w}^{*} S_{s w}^{*} + μ S_{s w}^{*}, \\ λ_{s w}^{*} S_{s w}^{*} & = (μ + γ) I_{s w}^{*}, \\ γ I_{s w}^{*} & = (μ + d) A_{s w}^{*}, \\ Λ_{c} & = λ_{c}^{*} S_{c}^{*} + μ S_{c}^{*}, \\ λ_{c} S_{c}^{*} & = (μ + γ) I_{c}^{*}, \\ γ I_{c}^{*} & = (μ + d) A_{c}^{*} \\ Λ_{f} & = λ_{f}^{*} S_{f}^{*} + μ S_{f}^{*}, \\ λ_{f} S_{f}^{*} & = (μ + γ) I_{f}^{*}, \\ γ I_{f}^{*} & = (μ + d) A_{f}^{*}, \\ Λ_{m} & = λ_{m}^{*} S_{m}^{*} + μ S_{m}^{*}, \\ λ_{m} S_{m}^{*} & = (μ + γ) I_{m}^{*}, \\ γ I_{m}^{*} & = (μ + d) A_{m}^{*} . \end{array} \end{array}

(11)

Where:

λ_{s w}^{*} = p_{M} (1 - π) \times \frac{Γ_{1} ψ_{s w}^{I_{c}} I_{c}^{*} + Γ_{2} ψ_{s w}^{A_{c}} A_{c}^{*}}{N_{c}},

λ_{c}^{*} = p_{F} [(1 - π) \times \frac{Γ_{1} ψ_{c}^{I_{s w}} I_{s w}^{*} + Γ_{2} ψ_{c}^{A_{s w}} A_{s w}^{*}}{N_{s w}} + \frac{Γ_{1} ψ_{c}^{I_{f}} I_{f}^{*} + Γ_{2} ψ_{c}^{A_{f}} A_{f}^{*}}{N_{f}}],

λ_{f}^{*} = p_{M} [\frac{Γ_{1} ψ_{f}^{I_{c}} I_{c}^{*} + Γ_{2} ψ_{f}^{A_{c}} A_{c}^{*}}{N_{c}} + \frac{Γ_{1} ψ_{f}^{I_{m}} I_{m}^{*} + Γ_{2} ψ_{f}^{A_{m}} A_{m}^{*}}{N_{m}}],

λ_{m}^{*} = p_{F} \times \frac{Γ_{1} ψ_{m}^{I_{f}} I_{f}^{*} + Γ_{2} ψ_{m}^{A_{f}} A_{f}^{*}}{N_{f}} .

Let us define the function $Ψ$ on $R_{+}^{*}$ by:

\begin{array}{r} Ψ (x) = x - 1 - l n (x) \end{array}

(12)

The function $Ψ$ is non-negative for all $x \in R_{+}^{*}$ . Let us consider the Lyapunov candidate function define by:

where:

\begin{array}{r} V_{s w} = & V_{S_{s w}} + V_{I_{s w}} + V_{A_{s w}}, \\ V_{c} = & V_{S_{c}} + V_{I_{c}} + V_{A_{c}}, \\ V_{f} = & V_{S_{f}} + V_{I_{f}} + V_{A_{f}}, \\ V_{m} = & V_{S_{m}} + V_{I_{m}} + V_{A_{m}} . \end{array}

and

V_{S_{s w}} = S_{s w}^{*} Ψ (\frac{S_{s w}}{S_{s w}^{*}}), V_{I_{s w}} = I_{s w}^{*} Ψ (\frac{I_{s w}}{I_{s w}^{*}}), V_{A_{s w}} = A_{s w}^{*} Ψ (\frac{A_{s w}}{A_{s w}^{*}}),

V_{S_{c}} = S_{c}^{*} Ψ (\frac{S_{c}}{S_{c}^{*}}), V_{I_{c}} = I_{c}^{*} Ψ (\frac{I_{c}}{I_{c}^{*}}), V_{A_{c}} = A_{c}^{*} Ψ (\frac{A_{c}}{A_{c}^{*}}),

V_{S_{f}} = S_{f}^{*} Ψ (\frac{S_{f}}{S_{f}^{*}}), V_{I_{f}} = I_{f}^{*} Ψ (\frac{I_{f}}{I_{f}^{*}}), V_{A_{f}} = A_{f}^{*} Ψ (\frac{A_{f}}{A_{f}^{*}})

V_{S_{m}} = S_{m}^{*} Ψ (\frac{S_{m}}{S_{m}^{*}}), V_{I_{m}} = I_{m}^{*} Ψ (\frac{I_{m}}{I_{m}^{*}}), V_{A_{m}} = A_{m}^{*} Ψ (\frac{A_{m}}{A_{m}^{*}}) .

Let us compute ${\dot{V}}_{S_{s w}}$

\begin{array}{r} {\dot{V}}_{S_{s w}} & = & (1 - \frac{S_{s w}^{*}}{S_{s w}}) {\dot{S}}_{s w} \\ = & (1 - \frac{S_{s w}^{*}}{S_{s w}}) [Λ_{s w} - \frac{S_{s w}}{N_{c}} p_{M} (1 - π) (Γ_{1} ψ_{s w}^{I_{c}} I_{c} + Γ_{2} ψ_{s w}^{A_{c}} A_{c} - μ S_{s w})] \end{array}

By using the first equation of (11), we have

\begin{array}{r} {\dot{V}}_{S_{s w}} & = & - μ \frac{{(S_{s w} - S_{s w}^{*})}^{2}}{S_{s w}} + \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{s w}^{*}}{N_{c}} [(1 - \frac{S_{s w}^{*}}{S_{s w}}) (1 - \frac{S_{s w} I_{s w}}{S_{s w}^{*} I_{s w}^{*}})] \\ + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{s w}^{*}}{N_{c}} [(1 - \frac{S_{s w}^{*}}{S_{s w}}) (1 - \frac{S_{s w} A_{s w}}{S_{s w}^{*} A_{s w}^{*}})] \\ {\dot{V}}_{S_{s w}} & = & - μ \frac{{(S_{s w} - S_{s w}^{*})}^{2}}{N_{c}} + \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{s w}^{*}}{N_{c}} (1 - \frac{S_{s w} I_{c}}{S_{s w}^{*} I_{c}^{*}} - \frac{S_{s w}^{*}}{S_{s w}} + \frac{I_{c}}{I^{*_{c}}}) \\ + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{s w}^{*}}{N_{c}} (1 - \frac{S_{s w} A_{c}}{S_{s w}^{*} A_{c}^{*}} - \frac{S_{s w}^{*}}{S_{s w}} + \frac{A_{c}}{A^{*_{c}}}) \\ = & - μ \frac{{(S_{s w} - S_{s w}^{*})}^{2}}{N_{c}} + \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{s w}^{*}}{N_{c}} (- \frac{S_{s w} I_{c}}{S_{s w}^{*} I_{c}^{*}} + 1 - l n (\frac{S_{s w} I_{c}}{S_{s w}^{*} I_{c}^{*}}) - \frac{S_{s w}^{*}}{S_{s w}} + 1 + l n (\frac{S_{s w}^{*}}{S_{s w}}) + \frac{I_{c}}{I^{*_{c}}} - 1 - l n (\frac{I_{c}}{I_{c}^{*}})) \\ + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{s w}^{*}}{N_{c}} (- \frac{S_{s w} A_{c}}{S_{s w}^{*} A_{c}^{*}} + 1 - l n (\frac{S_{s w} A_{c}}{S_{s w}^{*} A_{c}^{*}}) - \frac{S_{s w}^{*}}{S_{s w}} + 1 + l n (\frac{S_{s w}^{*}}{S_{s w}}) + \frac{A_{c}}{A^{*_{c}}} - 1 - l n (\frac{A_{c}}{A^{*_{c}}})) \end{array}

\begin{array}{r} {\dot{V}}_{S_{s w}} & = & - μ \frac{{(S_{s w} - S_{s w}^{*})}^{2}}{N_{c}} + \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{s w}^{*}}{N_{c}} [Ψ (\frac{S_{s w} I_{c}}{S_{s w}^{*} I_{c}^{*}}) - Ψ (\frac{S_{s w}}{S_{s w}^{*}}) + Ψ (\frac{I_{s w}}{I_{s w}^{*}})] \begin{array}{r} + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{s w}^{*}}{N_{c}} [Ψ (\frac{S_{s w} A_{c}}{S_{s w}^{*} A_{c}^{*}}) - Ψ (\frac{S_{s w}}{S_{s w}^{*}}) + Ψ (\frac{A_{s w}}{A_{s w}^{*}})] . \end{array} \end{array}

(13)

Let us compute ${\dot{V}}_{I_{s w}}$

\begin{array}{r} {\dot{V}}_{I_{s w}} & = & (1 - \frac{I_{s w}^{*}}{I_{s w}}) {\dot{I}}_{s w} \\ = & (1 - \frac{I_{s w}^{*}}{I_{s w}}) (\frac{S_{s w}}{N_{c}} (p_{M} (1 - π) (Γ_{1} ψ_{s w}^{I_{c}} I_{c} + Γ_{2} ψ_{s w}^{A_{c}} A_{c}) - (μ + γ) I_{s w}) \end{array}

From the second equation of system (11), we get:

\begin{array}{r} {\dot{V}}_{I_{s w}} & = & \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{c}^{*}}{N_{c}} (1 - \frac{I_{s w}^{*}}{I_{s w}}) (\frac{S_{s w} I_{c}}{S_{s w}^{*} I_{c}^{*}} - \frac{I_{s w}}{I_{s w}^{*}}) + \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{c}^{*}}{N_{c}} (1 - \frac{A_{s w}^{*}}{A_{s w}}) (\frac{S_{s w} A_{c}}{S_{s w}^{*} A_{c}^{*}} - \frac{A_{s w}}{A_{s w}^{*}}) \\ = & \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{c}^{*}}{N_{c}} (\frac{S_{s w} I_{s w}}{S_{s w}^{*} I_{s w}^{*}} - 1 - l n (\frac{S_{s w} I_{s w}}{S_{s w}^{*} I_{s w}^{*}}) - \frac{I_{s w}}{I_{s w}^{*}} + 1 + l n (\frac{I_{s w}}{I_{s w}^{*}}) - \frac{I_{s w}^{*} S_{s w} I_{c}}{I_{s w} S_{s w}^{*} I_{c}} 1 + l n (\frac{I_{s w}^{*} S_{s w} I_{c}}{I_{s w} S_{s w}^{*} I_{c}})) \\ + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{c}^{*}}{N_{c}} (\frac{S_{s w} A_{s w}}{S_{s w}^{*} A_{s w}^{*}} - 1 - l n (\frac{S_{s w} A_{s w}}{S_{s w}^{*} A_{s w}^{*}}) - \frac{A_{s w}}{A_{s w}^{*}} + 1 + l n (\frac{A_{s w}}{A_{s w}^{*}}) - \frac{A_{s w}^{*} S_{s w} A_{c}}{A_{s w} S_{s w}^{*} A_{c}} + 1 + l n (\frac{A_{s w}^{*} S_{s w} A_{c}}{A_{s w} S_{s w}^{*} A_{c}})) \end{array}

\begin{array}{r} {\dot{V}}_{I_{s w}} & = & \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{c}^{*}}{N_{c}} [Ψ (\frac{S_{s w} I_{c}}{S_{s w}^{*} I_{c}^{*}}) - Ψ (\frac{I_{s w}}{I_{s w}^{*}}) + Ψ (\frac{I_{s w}^{*} S_{s w} I_{c}}{I_{s w} S_{s w}^{*} I_{c}^{*}})] \\ + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{c}^{*}}{N_{c}} [Ψ (\frac{S_{s w} A_{c}}{S_{s w}^{*} A_{c}^{*}}) - Ψ (\frac{A_{s w}}{A_{s w}^{*}}) + Ψ (\frac{A_{s w}^{*} S_{s w} A_{c}}{A_{s w} S_{s w}^{*} A_{c}^{*}})] \end{array}

(14)

By adding the equation (13) and (14) we obtain

\begin{array}{r} {\dot{V}}_{I_{s w}} + {\dot{V}}_{I_{s w}} & = & - μ \frac{{(S_{s w} - S_{s w}^{*})}^{2}}{S_{s w}} + \frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{c}^{*}}{N_{c}} [- Ψ (\frac{S_{s w}}{S_{s w}^{*}}) + Ψ (\frac{I_{c}}{I_{c}^{*}}) - Ψ (\frac{I_{s w}}{I_{s w}^{*}}) + Ψ (\frac{I_{s w}^{*} S_{s w} I_{c}}{I_{s w} S_{s w}^{*} I_{c}^{*}})] \\ + & \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{c}^{*}}{N_{c}} [- Ψ (\frac{S_{c}}{S_{c}^{*}}) + Ψ (\frac{A_{c}}{A_{c}^{*}} + Ψ (\frac{I_{s w}^{*} S_{s w} A_{c}}{I_{s w} S_{s w}^{*} A_{c}^{*}})] \end{array}

(15)

Let us compute ${\dot{V}}_{A_{s w}}$

\begin{array}{r} {\dot{V}}_{A_{s w}} & = & (1 - \frac{A_{s w}^{*}}{A_{s w}}) {\dot{A}}_{s w} \\ = & (1 - \frac{A_{s w}^{*}}{A_{s w}}) (γ I_{s w} - (μ + d) A_{s w}) \\ = & γ I_{s w}^{*} (\frac{I_{s w}}{I_{s w}^{*}} - \frac{A_{s w}}{A_{s w}^{*}} - \frac{A_{s w}^{*} I_{s w}}{A_{s w} I_{s w}^{*}} + 1) \end{array}

\begin{array}{r} {\dot{V}}_{A_{s w}} & = & γ I_{s w}^{*} [\frac{I_{s w}}{I_{s w}^{*}} - \frac{A_{s w}}{A_{s w}^{*}} - \frac{A_{s w}^{*} I_{s w}}{A_{s w} I_{s w}^{*}}] \end{array}

(16)

Let $κ_{s w} = m a x {\frac{p_{M} (1 - π) Γ_{1} ψ_{s w}^{I_{c}} S_{s w}^{*} I_{c}^{*}}{N_{c}}; \frac{p_{M} (1 - π) Γ_{2} ψ_{s w}^{A_{c}} S_{s w}^{*} A_{c}^{*}}{N_{c}}; γ I_{s w}^{*}}$ .

From the relation (14), (15), (16) and $κ_{s w}$ , we get

\begin{array}{r} {\dot{V}}_{s w} & \leq & - μ \frac{{(S_{s w} - S_{s w}^{*})}^{2}}{S_{s w}} + κ_{s w} [- 2 Ψ (\frac{S_{s w}}{S_{s w}^{*}}) + Ψ (\frac{I_{c}}{I_{c}^{*}}) + Ψ (\frac{A_{c}}{A_{c}^{*}}) - Ψ (\frac{A_{s w}}{A_{s w}^{*}}) - Ψ (\frac{A_{s w}^{*} I_{s w}}{A_{s w} I_{s w}^{*}}) \\ + & Ψ (\frac{I_{s w}^{*} A_{s w} I_{c}}{I_{s w} S_{s w}^{*} I_{c}^{*}}) + Ψ (\frac{I_{s w}^{*} A_{s w} A_{c}}{I_{s w} S_{s w}^{*} A_{c}^{*}})] \end{array}

(17)

Consequently, by using the assumption $(H)$ , we obtain ${\dot{V}}_{s w} \leq 0$ . We are taking the same approach to show successively that ${\dot{V}}_{c} \leq 0$ , ${\dot{V}}_{m} \leq 0$ and ${\dot{V}}_{f} \leq 0$ . Hence $\dot{V} \leq 0$ . Thanks to LaSalle’s invariance principle, the solution $ε^{*}$ is said to be globally asymptomatically stable. This completes the proof.

Application with numerical simulations4.1. Parameters and initial data estimation

A population of 1,000,000 people is considered. According to the gender distribution in Burkina Faso [ 9], this makes 517,000 (51.7%) women and 483,000 men. Then, the proportion of FSWs among the female adults is considered to 0.85% [ 11]. Then the clients are assumed to represent 5% of the adult male population [ 25]. The number of infected men in both groups is calculated considering the prevalence of 0.5% [ 10]. The number of infected women is calculated considering the prevalence of HIV (5.4% among FSWs [ 11] and 0.8% in general population [ 10]).Table 2 provides the initial values for each compartment.

The natural mortality rate was estimated at µ = 0.0161 per year, based on life expectancy in Burkina Faso [ 9]. The other parameters were collected from literature or, when not available, just assumed. National estimates of paid sex between men are available from the 2010 Demographic and Health Survey [ 25]. However, this information does not appear in the 2021 version of the survey [ 26], which justifies the use of 2010 data for this analysis. Constant recruitment rates were arbitrarily chosen, considered at 2.94% of the population [ 9] at the initiation of the simulation.Table 3 gives the summary of parameters used for simulation.

4.2. Assessment of the impact of intervention

The intervention considered is a transmission blocking action such as condom use or pre-exposition prophylaxis. As these infections are rarely effective at 100%, several levels of effectiveness are considered. The parameter π is used to represent this effectiveness and five different levels are considered (Table 4).

Each scenario was simulated over 20 years and the impact of HIV prevention in the intervention population on the control population was estimated in terms of the annual number of incident HIV infections.

Table 2

Table 2. Initial condition ranges and baseline values used in the simulation

Compartment	Baseline	Compartment	Baseline
S_c(0)	24,029	S_sw(0)	4,162
I_c(0)	121	I_sw(0)	238
A_c(0)	1	A_sw(0)	1
S_m(0)	456,556	S_f (0)	508,499
I_m(0)	2,294	I_f (0)	4,101
A_m(0)	1	A_f (0)	1

Table 3

Table 3. Parameter ranges and baseline values per unit year.

Parameter	Min	Max	Baseline	Source	Parameter	Min	Max	Baseline	Source
death	-	-	0.01512	[25]	-	-	-	-	-

Table 4

Table 4. Levels of effectiveness considered for the interventions.

Values of	Comments
0.00	No intervention at all
0.20	Intervention with low efficacy
0.50	Intervention with medium efficacy
0.80	Intervention with high efficacy
1.00	Intervention with full efficacy

4.3. Numerical results

The populations simulated are presented inFigure figure 2 prior to any intervention. The red line represents the situation in absence of any intervention while the black line represents an ideal situation where the transmission is completely stopped between FSWs and their clients. The violet, blue and green lines will represent the different levels of efficacy of the intervention program. In the absence of any intervention, the high infection rate among the FSWs brings to a very fast decrease in the population of susceptible. This population is divided into two in one year’s time and continues at this pace exponentially. The population of clients free of infection is also decreasing. Even though that decrease is not exponential, it is also fast with an increasing percentage of cases yearly. The general population however presents a linear trend.

While looking at the intervention, its effect on the target population is clearly observed. The number of people in the susceptible compartment is increasing as the intervention is effective. In the case of the ideal intervention with 100% efficacy, infections and transmissions are completely blocked between FSWs and their clients. Then, the susceptible cases in this population are no longer reducing but it is even increasing due to the natural population growth.

Figure 2

Evolution of the number of individuals free of infection by category for the different levels of intervention.

The effectiveness of the intervention can be clearly observed inFigure figure 3, presenting the evolution of the people in the infectious class. In absence of intervention, the epidemic curve of the FSWs presents the classical shape where a high increase is observed, until a peak and then a slow decrease occurs due to the drastic reduction of the susceptible population. This leads to the extinction of the population if nothing is done. Then, high impacts are observed among the FSWs when we consider the graphic. The more efficient the intervention is, the lower the curve is. Finally, looking at the fully efficient curve, the infectious population just reduces to disappearing when the infectious cases before the intervention will all move to AIDS state. The effect of the intervention is also observed among clients as the FSWs are considered their principal infecting agents. Once again, the intervention will reduce the growth speed and when it is fully efficient, it can even lead to a reduction in the number of infectious diseases.

The principal aspect which is generally hard to observe is the effect of the intervention in the general population. ObservingFigure figure 3, the epidemic curves are superimposed, showing a non-observable effect. However, the curves present a clearly decreasing trend showing that the infection is well controlled in this general population.

Figure 3

Evolution of the number of cases in infectious state by category for the different level of intervention

Figure 4 presents the difference in the number of infected persons between intervention scenarios and no intervention at all. This provides a clear idea on the number of infections avoided every year. At this moment, it is possible to observe the effect of the intervention in the general population. For example, for the simulation performed, after 10 years, in the female population, the equivalent of 0.7 cases are saved for the low impact simulation. More than 3 infections are saved for the high impact intervention and for the full efficient program. In the full impact intervention, 35 infections are avoided after 20 years. The impact of the intervention is on the other hand very low in among general men population. This is understandable as the men do not have direct interaction with the intervention target group as the women have interaction with clients.

Conclusion

This study presented a compartmental model of HIV transmission between high-risk populations (FSWs and clients) and low-risk populations (men and women in the general population). It considered prevention action among the high-risk population. The model was calibrated according to epidemiological and demographic parameters reflecting the context of Burkina Faso and only sexual HIV transmission was considered. The mathematical analysis of the model was carried out before proceeding with the numerical simulation. This involved checking the basic properties of the model, namely the positivity of solutions, boundedness, and calculating the basic reproduction number ( $R_{0}$ ) associated with the model. Using Varga’s theorem, it was shown that the disease-free equilibrium ε₀ was globally asymptotically stable when $R_{0} < 1$ . The Lyapunov function technique was used to obtain the global stability of the endemic equilibrium ε^∗, when $R_{0}$ > 1. Finally, to simulate the impact of HIV prevention among FSWs and their clients on HIV transmission in the general population, different simulation scenarios were tested by varying the efficacy of the intervention. The results of the model suggest that these parameters play a major role in reducing the number of people infected with HIV in both sub-populations. The evolution of the infection in all sub populations, compared with a status quo situation where no intervention is performed can be observed. Then, the impact of the intervention can be estimated as risk ratio or incidence rate ratio.

Figure 4

Evolution of the number of infections prevented by the intervention by category for the different level of intervention.

This work should not be considered as an actual estimation of a given intervention targeting FSWs in Burkina Faso. The intervention proposed here is totally hypothetical and arbitrary values were taken for the efficiency level. Additionally, the input parameters are still questionable, and one could prefer another source. This is to be considered as proof of concept. In fact, it is possible, with applied mathematical tools, to go beyond the estimation of the direct effect of an intervention and find the long-term effect, concerning everyone. This kind of estimation will bring more arguments to promote HIV control programs on key populations. These persons are often discriminated against or even condemned by society. It is therefore often difficult to justify the collective efforts needed to reduce infection at their level. However, by applying the approach developed here, it is possible to see the effect of these efforts on any member of society. In other words, we demonstrate that it is in every individual’s interest to protect key populations from transmission. This approach can be generalized to any context, any country, provided the basic data are availables. Above all, it can be extended to other key populations (MSM, transgender people, prisoners, etc.), provided that the disease transmission pattern is redefined.

However, it is important to stress that the model presented in this work is a highly simplified caricature of a complex interaction between at-risk populations and therefore has certain obvious limitations. For example, the model does not take into account the abandonment of risky behaviours by FSWs and their clients. In another hand, no prevention action is considered for general population. The experimental data needed to verify the model is also scarce, particularly data on FSW clients. There are several ways of getting around these limitations. First, it is necessary to link the model to the data for a clearer estimation of the parameter values, a Bayesian approach for example. Nevertheless, despite the limitations highlighted, the results of the model have a significant impact on the dynamics of HIV in at-risk populations.

Acknowledgements: The author thank the editor and referees for their careful reading and valuable comments which led to a significant improvement of the original manuscript.

Conflict of interest: The authors declare no conflict of interest.

References [1]

ONUSIDA. StratÉgie mondiale de lutte contre le sida, 2021-2026 mettre fin aux inÉgalitÉs. mettre fin au sida., 2021. https://www.unaids.org/sites/default/files/media_asset/ global-AIDS-strategy-2021-2026_fr.pdf.

[2]

UNAIDS. Unaids data 2023, 2024. https://www.unaids.org/sites/default/files/media_asset/ data-book-2023_en.pdf.

[3]

ONUSIDA. Regional factsheet western and central africa, 2023. https://thepath.unaids.org/ wp-content/themes/unaids2023/assets/files/regional_fs_western_central_africa.pdf.

[4]

Joint United Nations Programme on HIV/AIDS (UNAIDS). Prevention gap report. Technical report, UNAIDS, September 2016.

[5]

EO Omondi, RW Mbogo, and LS Luboobi. A mathematical model of hiv transmission between commercial sex workers and injection drug users. Research in Mathematics, 9(1):1-15, 2022.

[6]

M Maheu-Giroux, J F Vesga, S Diabaté, M Alary, S Baral, D Diouf, K Abo, and M-C Boily. Changing dynamics of hiv transmission in côte d'ivoire: modeling who acquired and transmitted infections and esti- mating the impact of past hiv interventions (1976-2015). Journal of acquired immune deficiency syndromes (1999), 75(5):517, 2017.

[7]

L Geidelberg, K M Mitchell, M Alary, A Mboup, L Béhanzin, F Guédou, N Geraldo, E Goma-Matsétsé, K Giguère, M Aza-Gnandji, et al. Mathematical model impact analysis of a real-life pre-exposure prophy- laxis and treatment-as-prevention study among female sex workers in cotonou, benin. Journal of Acquired Immune Deficiency Syndromes (1999), 86(2): e28, 2021.

[8]

A Low, N Nagot, I Konate, N Meda, M Segondy, P Van de Perre, P Mayaud, and P Vickerman. Potential impact of existing interventions and of antiretroviral use in female sex workers on transmission of hiv in burkina faso: a modeling study. JAIDS Journal of Acquired Immune Deficiency Syndromes, 68: S180-S188, 2015.

[9]

INSD. Cinquième recensement général de la population et de l'habitation du burkina faso. Technical report, Institut National de la Statistiqu et de la Démographie, 2019.

[10]

ONSP. Profil du burkina faso sur le vih-ist. Technical report, Observatoire National de la Santé de la Population, 2023. http://cns.bf/IMG/pdf/mshp_profil_du_burkina_sur_le_vih.pdf.

[11]

CNLS. Enquête biocomportementale auprès des populations clés au burkina faso. Technical report, Secr'étariat Permanent du Conseil National de Lutte contre le SIDA et les Infections Sexuellement Trans- missibles (SP-CNLS/IST), 2017.

[12]

S Saravanakumar, A Eswari, L Rajendran, and Marwan Abukhaled. A mathematical model of risk fac- tors in hiv/aids transmission dynamics: observational study of female sexual network in india. Applied Mathematics & Information Sciences, 14(6):967-976, 2020.

[13]

O Koutou, B Sangaré, et al. Mathematical analysis of the impact of the media coverage in mitigating the outbreak of covid-19. Mathematics and Computers in Simulation, 205:600-618, 2023.

[14]

R Bellman, K L Cooke, et al. Diflerential-diflerence equations. Academic press, 1963.

[15]

O Koutou, B Traoré, and B Sangaré. Mathematical modeling of malaria transmission global dynamics: taking into account the immature stages of the vectors. Advances in Diflerence Equations, 2018:1-34, 2018.

[16]

O Diekmann, J Andre Peter Heesterbeek, and J AJ Metz. On the definition and the computation of the basic reproduction ratio r 0 in models for infectious diseases in heterogeneous populations. Journal of mathematical biology, 28:365-382, 1990.

[17]

P Van den Driessche and J Watmough. Reproduction numbers and sub-threshold endemic 628 equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180(1-2):629, 2002.

[18]

Richard S Varga. Matrix iterative analysis: Englewood cliffs. XIV-NY: Prentice-Hall, 1962.

[19]

H Ouedraogo, D Ouedraogo, I Ibrango, and A Guiro. A study of stability of seihr model of infectious disease transmission. Nonautonomous Dynamical Systems, 8(1):307-327, 2021.

[20]

A Diabate, B Sangaré, and O Koutou. Optimal control analysis of a covid-19 and tuberculosis (tb) co- infection model with an imperfect vaccine for covid-19. SeMA Journal, 81, 06 2023.

[21]

A Kaboré, B Sangaré, and B Traoré. Mathematical model of mosquito population dynamics with constants and periodic releases of wolbachia-infected males. Applied Mathematics in Science and Engineering, 32, 07 2024.

[22]

B Sangaré. Contribution à la modélisation mathématique en épidémiologieet en écologie et à la résolution numérique des EDP évolutives. Habilitation à Diriger des Recherches, Université de Franche-Comté, June 2023.

[23]

V Lakshmikantham, S Leela, and A A Martynyuk. Stability analysis of nonlinear systems. Springer, 1989.

[24]

JP LaSalle. The stability of dynamical systems, society for industrial and applied mathematics, philadel- phia, pa., 1976. In With an appendix:"Limiting equations and stability of nonautonomous ordinary difler- ential equations" by Z. Artstein, Regional Conference Series in Applied Mathematics, 1976.

[25]

INSD. Enquête démographique et de santé et à indicateurs multiples. Technical report, Institut Nationale de la Statistiqu et de la Démographie, 2012.

[26]

INSD. Enquête démographique et de santé et à indicateurs multiples. Technical report, Institut Nationale de la Statistiqu et de la Démographie, 2022. https://dhsprogram.com/pubs/pdf/FR378/FR378.pdf.

[27]

Protéger les Canadiens contre les maladies. Risque de transmission du vih: Sommaire des données scien- tifiques, 2012.

[1]

[2]

UNAIDS. Unaids data 2023, 2024. https://www.unaids.org/sites/default/files/media_asset/ data-book-2023_en.pdf.

[3]

ONUSIDA. Regional factsheet western and central africa, 2023. https://thepath.unaids.org/ wp-content/themes/unaids2023/assets/files/regional_fs_western_central_africa.pdf.

[4]

Joint United Nations Programme on HIV/AIDS (UNAIDS). Prevention gap report. Technical report, UNAIDS, September 2016.

[5]

EO Omondi, RW Mbogo, and LS Luboobi. A mathematical model of hiv transmission between commercial sex workers and injection drug users. Research in Mathematics, 9(1):1-15, 2022.

[6]

[7]

[8]

[9]

INSD. Cinquième recensement général de la population et de l'habitation du burkina faso. Technical report, Institut National de la Statistiqu et de la Démographie, 2019.

[10]

ONSP. Profil du burkina faso sur le vih-ist. Technical report, Observatoire National de la Santé de la Population, 2023. http://cns.bf/IMG/pdf/mshp_profil_du_burkina_sur_le_vih.pdf.

[11]

[12]

[13]

O Koutou, B Sangaré, et al. Mathematical analysis of the impact of the media coverage in mitigating the outbreak of covid-19. Mathematics and Computers in Simulation, 205:600-618, 2023.

[14]

R Bellman, K L Cooke, et al. Diflerential-diflerence equations. Academic press, 1963.

[15]

[16]

[17]

P Van den Driessche and J Watmough. Reproduction numbers and sub-threshold endemic 628 equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180(1-2):629, 2002.

[18]

Richard S Varga. Matrix iterative analysis: Englewood cliffs. XIV-NY: Prentice-Hall, 1962.

[19]

H Ouedraogo, D Ouedraogo, I Ibrango, and A Guiro. A study of stability of seihr model of infectious disease transmission. Nonautonomous Dynamical Systems, 8(1):307-327, 2021.

[20]

A Diabate, B Sangaré, and O Koutou. Optimal control analysis of a covid-19 and tuberculosis (tb) co- infection model with an imperfect vaccine for covid-19. SeMA Journal, 81, 06 2023.

[21]

[22]

[23]

V Lakshmikantham, S Leela, and A A Martynyuk. Stability analysis of nonlinear systems. Springer, 1989.

[24]

[25]

INSD. Enquête démographique et de santé et à indicateurs multiples. Technical report, Institut Nationale de la Statistiqu et de la Démographie, 2012.

[26]

[27]

Protéger les Canadiens contre les maladies. Risque de transmission du vih: Sommaire des données scien- tifiques, 2012.

[1]

[2]

UNAIDS. Unaids data 2023, 2024. https://www.unaids.org/sites/default/files/media_asset/ data-book-2023_en.pdf.

[3]

ONUSIDA. Regional factsheet western and central africa, 2023. https://thepath.unaids.org/ wp-content/themes/unaids2023/assets/files/regional_fs_western_central_africa.pdf.

[4]

Joint United Nations Programme on HIV/AIDS (UNAIDS). Prevention gap report. Technical report, UNAIDS, September 2016.

[5]

EO Omondi, RW Mbogo, and LS Luboobi. A mathematical model of hiv transmission between commercial sex workers and injection drug users. Research in Mathematics, 9(1):1-15, 2022.

[6]

[7]

[8]

[9]

INSD. Cinquième recensement général de la population et de l'habitation du burkina faso. Technical report, Institut National de la Statistiqu et de la Démographie, 2019.

[10]

ONSP. Profil du burkina faso sur le vih-ist. Technical report, Observatoire National de la Santé de la Population, 2023. http://cns.bf/IMG/pdf/mshp_profil_du_burkina_sur_le_vih.pdf.

[11]

[12]

[13]

O Koutou, B Sangaré, et al. Mathematical analysis of the impact of the media coverage in mitigating the outbreak of covid-19. Mathematics and Computers in Simulation, 205:600-618, 2023.

[14]

R Bellman, K L Cooke, et al. Diflerential-diflerence equations. Academic press, 1963.

[15]

[16]

[17]

P Van den Driessche and J Watmough. Reproduction numbers and sub-threshold endemic 628 equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180(1-2):629, 2002.

[18]

Richard S Varga. Matrix iterative analysis: Englewood cliffs. XIV-NY: Prentice-Hall, 1962.

[19]

H Ouedraogo, D Ouedraogo, I Ibrango, and A Guiro. A study of stability of seihr model of infectious disease transmission. Nonautonomous Dynamical Systems, 8(1):307-327, 2021.

[20]

A Diabate, B Sangaré, and O Koutou. Optimal control analysis of a covid-19 and tuberculosis (tb) co- infection model with an imperfect vaccine for covid-19. SeMA Journal, 81, 06 2023.

[21]

[22]

[23]

V Lakshmikantham, S Leela, and A A Martynyuk. Stability analysis of nonlinear systems. Springer, 1989.

[24]

[25]

INSD. Enquête démographique et de santé et à indicateurs multiples. Technical report, Institut Nationale de la Statistiqu et de la Démographie, 2012.

[26]

[27]

Protéger les Canadiens contre les maladies. Risque de transmission du vih: Sommaire des données scien- tifiques, 2012.

[1]

[2]

UNAIDS. Unaids data 2023, 2024. https://www.unaids.org/sites/default/files/media_asset/ data-book-2023_en.pdf.

[3]

ONUSIDA. Regional factsheet western and central africa, 2023. https://thepath.unaids.org/ wp-content/themes/unaids2023/assets/files/regional_fs_western_central_africa.pdf.

[4]

Joint United Nations Programme on HIV/AIDS (UNAIDS). Prevention gap report. Technical report, UNAIDS, September 2016.

[5]

EO Omondi, RW Mbogo, and LS Luboobi. A mathematical model of hiv transmission between commercial sex workers and injection drug users. Research in Mathematics, 9(1):1-15, 2022.

[6]

[7]

[8]

[9]

INSD. Cinquième recensement général de la population et de l'habitation du burkina faso. Technical report, Institut National de la Statistiqu et de la Démographie, 2019.

[10]

ONSP. Profil du burkina faso sur le vih-ist. Technical report, Observatoire National de la Santé de la Population, 2023. http://cns.bf/IMG/pdf/mshp_profil_du_burkina_sur_le_vih.pdf.

[11]

[12]

[13]

O Koutou, B Sangaré, et al. Mathematical analysis of the impact of the media coverage in mitigating the outbreak of covid-19. Mathematics and Computers in Simulation, 205:600-618, 2023.

[14]

R Bellman, K L Cooke, et al. Diflerential-diflerence equations. Academic press, 1963.

[15]

[16]

[17]

P Van den Driessche and J Watmough. Reproduction numbers and sub-threshold endemic 628 equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180(1-2):629, 2002.

[18]

Richard S Varga. Matrix iterative analysis: Englewood cliffs. XIV-NY: Prentice-Hall, 1962.

[19]

H Ouedraogo, D Ouedraogo, I Ibrango, and A Guiro. A study of stability of seihr model of infectious disease transmission. Nonautonomous Dynamical Systems, 8(1):307-327, 2021.

[20]

A Diabate, B Sangaré, and O Koutou. Optimal control analysis of a covid-19 and tuberculosis (tb) co- infection model with an imperfect vaccine for covid-19. SeMA Journal, 81, 06 2023.

[21]

[22]

[23]

V Lakshmikantham, S Leela, and A A Martynyuk. Stability analysis of nonlinear systems. Springer, 1989.

[24]

[25]

INSD. Enquête démographique et de santé et à indicateurs multiples. Technical report, Institut Nationale de la Statistiqu et de la Démographie, 2012.

[26]

[27]

Protéger les Canadiens contre les maladies. Risque de transmission du vih: Sommaire des données scien- tifiques, 2012.