Analysing transmission dynamics of HIV/AIDS with optimal control strategy and its controlled state

HIV is a virus that weakens a person's immune system. HIV has three stages, and AIDS is the most severe stage of HIV (Stage 3). People with HIV should take medicine (called ART) recommended by WHO as soon as possible to reduce the amount of virus in the body. In this paper, we formulate a mathematical model for HIV/AIDS with a new approach by focusing on two groups of infectious individuals, HIV and AIDS. We also introduce a controlled class (treated patients and being monitored), and people in this class can spread the disease. We further investigate the essential dynamics of the model through an equilibrium analysis. Optimal control theory is applied to explore effective treatment strategies by combining two control measures: standard antiretroviral therapy and AIDS treatments. Numerical simulation results show the effects of the two time-dependent controls, and they can be used as guidelines for public health interventions.


Introduction
HIV (human immunodeficiency virus) is a virus that attacks the body's immune system. There are three stages of HIV: acute HIV infection, chronic HIV infection, and acquired immunodeficiency syndrome (AIDS). The last stage is the most severe phase of HIV infection. With proper treatments, some HIV-infected people may never move to the last stage. HIV remains a serious public health threat worldwide, especially in some developing countries and more than 54 % of 37 million infected people are aware of their infection [36]. In 2020, about 680,000 people died from HIV-related causes. According to CDC, having anal sex, sharing needles, syringes, or other drug-injection equipment are the main causes to pass HIV from one person to another [34]. Getting HIV from the workplace, pre-chewed food, oral sex, or biting are rare. HIV can be successfully controlled by the provision of personal protection such as using condoms, and people with HIV may never move into stage 3 (AIDS) if they take HIV medicine as prescribed. There is no cure for HIV infection, however, there are ways to control it. Staying healthy and taking HIV medicine as prescribed can make the viral load very low. The well-known HIV treatment is called antiretroviral therapy (ART), and it is recommended by WHO for those who are living with HIV. It consists of the combination of antiretroviral (ARV) drugs to maximally suppress HIV and stop the progression of HIV disease. Also, the reduced number of infections and deaths by using ART is significant. WHO suggested that ART should be deployed to infected individuals with HIV as soon as possible after diagnosis. By June 2021, WHO reported that 187 countries had adopted the recommendation. To better understanding the complex transmission dynamics of HIV/AIDS, several mathematical models have been proposed and analysed in recent years [2,7,8,11,15,24,26,27,31].
In 2009, Mukandavire et al. proposed a model to study the effects of public health education campaigns on the spread of HIV/AIDS as a single-strategy approach in HIV prevention [17]. They studied the effect of an education program for preventing the spread of HIV by comparing the basic reproductive number for the model with education measure to the basic reproduction number (R 0 ) for HIV/AIDS in the absence of any intervention to assess the possible community benefits of public health education campaigns.
In 2010, Mohamad Shirazian et al. formulated an optimal control mathematical model for HIV infection, and they also investigated the total cost of treatments [31]. The numerical results were implemented by using an AVK method. The paper mainly proposed an estimation method to evaluate parameters in the model, and the results are based on several assumptions. Also, the model focused on the number of CD4 T cells instead of infected people.
Later in 2011, Mukandavire et al. formulated an SIA model to study the HIV/AIDS epidemic trends in South Africa [19]. The results show that the model can predict the spread of the disease very well compared with other research outcomes on the HIV/AIDS epidemic in South Africa. In the same year, Mushayabasa et al. proposed a deterministic model to investigate the effects of socioeconomic status on the transmission dynamics of HIV/AIDS [18]. They studied two groups of people. The first group (risk group) was considered poor people. In contrast, rich people were classified as a low-risk group. Results suggested that if individuals referred to as poor involve risky sexual behaviour, HIV/AIDS persisted in the community. In the same year, Yan Wang et al. introduced a delay-dependent model with HIV drug resistance during therapy [32]. They studied the evolution of drug resistance by considering an HIV viral dynamical model with a general form of target cell density, drug resistance, and intracellular delay incorporating antiretroviral therapy. They established sufficient conditions for global asymptotic stability of the disease-free and drug-resistant steady states. The basic reproduction numbers for each case, drug-sensitive and drug-resistant strains were calculated. In addition, the study concludes that sustained oscillations from numerical simulations mainly depend on the target cell density rather than the intracellular delay.
Defang Liu et al. proposed a novel time-delayed HIV/AIDS model in 2013 [7]. They investigated the effect of vaccination and antiretroviral therapy on the model. Also, they introduced three states of AIDS according to the developing progress of infection. The results show that the therapy for HIV-infected-aware individuals could reduce the number of infective people. In the same year, Okosun et al. purposed two optimal control models [20,21]. The first one focused on the analysis of the recruitment and industrial human resources management for optimal productivity in the presence of HIV/AIDS epidemic. The effects of prevention, education, and HAART treatment on the spread of HIV/AIDS were investigated. The result shows that the effective education of employees and adherence to preventive measures against infection is the most cost-effective strategy. The second one they studied about the impact of optimal control on the treatment of HIV/AIDS and screening of unaware infectives. The study shows that the successful screening of unaware infectives has a significant impact in reducing the endemicity of HIV/AIDS.
In 2014, Pablo S. Rivadeneira et al. presented a model of HIV dynamics with antiretroviral therapy (ART) [27]. They provided ways to find treatment strategies by computing the drug treatment directly in amounts of drugs, which is very helpful to approximate the costs of public health interventions. In the same year, Takaidaza et al. purposed an optimal control of HIV/AIDS model that aimed to investigate the disease transmission in the community with substance abuse problem [33]. The study explored the implementation of optimal control strategies involving treatment of substance-abusing susceptibles, counselling, and prevention to fight the spread of HIV.
Mojaver et al. formulated an HIV infection model that studied the cell-free virus diffusion and direct cell-to-cell transmission in 2015 [16]. The study was based on some important biological meanings of HIV infection among cells. They also included antiretroviral therapy in their study. Later, Seidu et al. formulated optimal control models of HIV/AIDS and HIV-Malaria in the workplace [28][29][30]. For HIV/AIDS model, the study aimed to reduce infections in the workplace by combining some controls: reducing the infection of susceptible individuals, treating infected individuals, changing behaviour, and reducing non-productivity at the workplace. The study shows that employing all control efforts is most effective in the fight against the disease. Similar to the HIV-Malaria model, implementing the two preventive measures, Malaria and HIV/AIDS prevention controls, is the most effective strategy. In 2017, Nigar Ali et al. [1] presented an optimal control strategy of the HIV-1 epidemic model for the recombinant virus. They implemented two control variables: drug therapy in blocking off the infection of new cells and drug therapy in decreasing virus. However, the study lacked a rigorous analysis of the system's stability. The study concludes that the viral load reduces after treatment, and the CD4T population increases. This result shows that the infected cells drop, and the healthy cells increase.
Later in 2018, Omondi et al. [22] formulated an HIV/AIDS model that studied the impact of testing, treatment, and control of HIV transmission in Kenya. The results show that the combination of controls could reduce the spread of the disease. However, the control rates were high for most of the epidemic period, and some developing countries could not afford this plan. The main focus of this study was to investigate the transmission dynamics between two infectious classes, lower plasma viral load, and high plasma viral load.
Kheiri et al. presented a fractional model for transmission of HIV/AIDS in 2018 [13]. They also investigated the effects of condom use and antiretroviral therapy (ART) controls [12]. The study represented the effects of using ART as an HIV control, and the results show that the numbers of diagnosed, undiagnosed HIV-infected population, and AIDS people are reduced. These results confirm that the ART treatment is an effective strategy.
In 2019, Omondi et al. [23] continued their work on the modelling of HIV infection in two heterosexual age groups in Kenya. The study focused on HIV controls to guide public interventions: condom distribution and antiretroviral therapy. They used a simple SIT model (Susceptible, Infectious, and ART) to understand the dynamics of HIV within and between two different age groups in Kenya. The results show that most HIV infections occur in the adult population.
From these works, we see that ART is an effective medicine to control the disease. Taking ART as prescribed can also reduce the risk that infected people spread HIV to others. There are two main interests in those studies. One studies cells, and another focuses on people. In this paper, we formulate an HIV-AIDS mathematical model with a new approach by focusing on two groups of infectious people: HIV and AIDS. We also introduce a controlled class (patients are treated and being monitored) into the model. The transmission rates among these states are carefully considered. The study also includes testing the robustness of the proposed model as a simple and effective tool for understanding HIV transmission dynamics and controls. Control measures using ART and medical treatment for AIDS patients are investigated. We begin our study by analysing an autonomous version of the model system with constant controls and then carry out a rigorous mathematical analysis of our model. Further, we consider a model system with time-dependent controls to examine the optimal implementation of the intervention measures by assessing their effects and costs.
In the remainder of this paper, we will first present the HIV-AIDS model with control measures and its controlled class of HIV patients. We will then conduct an equilibrium analysis for the epidemic and endemic dynamics of the system when the controls are constants. Then, we will turn to a time-dependent control system and perform an optimal control study for the HIV-AIDS model. Finally, we round up the paper by conclusion and discussion.

HIV-AIDS model
In this section, we derive a mathematical model with antiretroviral therapy and treatment for AIDS patients. The model classifies the human population, denoted by N, into the susceptibles (S), the infected or infectives with HIV (I), the AIDS patients (A), and the controlled (C). We use the letter C for the controlled class to indicate that these people have had HIV with no symptoms, but the virus still exists in their blood. Thus, this class is considered as an almost recovered class from HIV. Also, people in this class are aware of having HIV, so they are more careful to prevent themselves from spreading the disease to others.
We assume that individuals are born and die at an average rate of μ. We also denote the disease-caused death rates for HIV and AIDS patients by α 1 and α 2 , respectively. Susceptible individuals acquire HIV through HIV individuals, AIDS patients, or controlled individuals, at rates β 1 , β 2 , and β 3 , respectively. We further assume that the infected individuals with HIV are treated with ART at a rate of φ 1 (t), and AIDS patients are treated medical treatments at a rate φ 2 (t). The treated patients are moved into the controlled class. Also, the infected with HIV class moves into AIDS class at a rate of σ . Our model thus takes the form below (see Figure 1 for its diagram): (1) In case there is no disease-related mortality, we have N = S + I + A + C the constant total population. In general, φ 1 (t) and φ 2 (t) are functions of t, representing non-uniform and timedependent controls. For the special case when the rates of all the two controls are positive constants, i.e. φ 1 (t) = φ 1 > 0 and φ 2 (t) = φ 2 > 0. The model (1) is reduced to an autonomous system. This allows us to conduct a careful equilibrium analysis, in Sections 3 and 4, on both the disease-free and endemic equilibria. It is clear to see that in this case, the feasible region of the system (1) is which is positively invariant for the flow of the system (1).

Disease-free equilibrium
With constant controls, the disease-free equilibrium (DFE) of the system (1) is given by To compute the basic reproductive number, R 0 , we use the well-known method of van den Driessche and Watmough [35]. Here the associated next-generation matrices are given by The basic reproductive number is then determined as the spectral radius of FV −1 , which yields Consequently, based on the work in the paper proposed by van den Driessche and Watmough [35], we immediately have the following result: Indeed, the DFE is also globally asymptotically stable when R 0 < 1. To study the global asymptotic stability of DFE, one common approach is to construct an appropriate Lyapunov function. We have found, however, that it is simpler to apply the following result introduced by Castillo-Chavez et al. [6].

Lemma 3.1: Consider a model system written in the form
where X 1 ∈ R m denotes (its components) the number of uninfected individuals and X 2 ∈ R n denotes (is components) the number of infected individuals including latent, infectious, etc; we denote the disease-free equilibrium of the system by X 0 = X * 1 . Also assume the conditions (H1) and (H2) below: is an M-matrix (the off diagonal elements of A are non-negative) and is the region where the model makes biological sense. Then the DFE X 0 = (X * 1 , 0) is globally asymptotically stable provided that R 0 < 1.

Theorem 3.2:
When R 0 < 1, the DFE of the system (1) with constant controls is globally asymptotically stable.

Proof:
We only need to show that the condition (H1) and (H2) hold when R 0 < 1. In our ODE system, X 1 = S, X 2 = (I, A, C), and X * 1 = N. We note that the system is linear and its solution can be easily found as We have dX 1 dt = F(X 1 , 0) = [μN − μS], is linear and its solution can be easily found as . Clearly, S(t) → N as t → ∞, regardless of the values of S(0). Thus X * 1 = N is globally asymptotically stable. Now, consider where A = ∂G ∂X 2 (X * 1 , 0). Specifically, where the off-diagonal elements are obviously non-negative. Substituting into (5) givesĜ(X 1 , (2) holds as well and we complete the proof.

Endemic equilibrium
When the disease is present in the population, I * = 0, there may be several critical points where I * = 0, which are the endemic equilibrium points (EEP) of the model. These points are denoted by * = (S * , I * , A * , C * ) which are determined from the model as follows: Its components must satisfy, We first show the following theorem.
Substitute of A * , C * and S * in I * , we have With some manipulations, we now have an equation for I * as follows: which can be written in the form H 1 I * 2 + H 2 I * = 0, where From the reproduction number Thus we conclude that I * > 0 and we complete the proof.

Local and global stabilities
We proceed to analyse the stability properties of the endemic equilibrium. First, we prove the following result regarding local stability.

Theorem 4.2:
The positive endemic equilibrium * is locally asymptotically stable.
Proof: The Jacobian of the system (1) at x = * is given by We choose 3rd row implies that Conformational the above equation, we obtain that giving the following characteristic polynomial, where To ensure that all roots of Equation (6) have negative real parts, the Routh-Hurwitz stability criterion will be used. We require that a 1 > 0, a 3 > 0anda 1 a 2 − a 3 > 0. Now, we note that We can arrange the equation dI * dt = 0 in the other form as follows: Thus We see that Equations (7) and (8) are greater than to zero, this implies that Next, consider that Let dI * dt = 0, we have Taking μh 2 × (9), we obtain Setting dA * dt = 0, we have σ I * − (μA * + φ 2 A * ) = 0. That is Substitute of Equation (11) in Equation (10), we have Next, let dS * dt = 0. Then Taking h 1 h 2 × (13), we obtain Thus Adding (12) + (14); Indeed, by setting dS * dt = 0 = dI * dt . We obtain Taking (15) + (16), we have μN − μS * − h 1 I * = 0. This implies that a 3 > 0, as required.
Next, we will show that a 1 a 2 > a 3 ⇔ a 1 a 2 − a 3 > 0. Consider Bringing the common factor (M − W + μ + h 1 + β 3 S * − β 1 S * ) out of the above equation, we obtain We will show that a 1 a 2 − a 3 > 0 by verifying each of components is greater than zero as follows: Clearly, M − W + μ + h 1 + β 3 S * − β 1 S * = a 1 > 0. Consider the last component of the equation which is always true. Next, we will show that Since a 1 > 0 and a 3 > 0, we can obtain a 3 h 2 and h 2 × a 1 , respectively as the following: Taking (18) + (19), we have We observe that Equations (17) and (20) are similar, so we take different terms to compare, implies that σ ( That is h 2 2 + a 2 > 0. Hence a 1 a 2 − a 3 > 0. We complete the proof.
Next, we will follow the geometric approach originally proposed by Li and Muldowney [9] to investigate the global asymptotic stability of the endemic equilibrium. To that end, we first present the following result based on the geometric approach.
In Equation (21), P is a matrix-valued function defined as where Q(X) is a ( n 2 ) × ( n 2 ) matrix-valued C 1 function in D, Q f is the derivative of Q (entrywise) along the direction of f , and J [2] is the second additive compound matrix of the Jacobian J(X) = Df (X). Meanwhile, m(P) is the Lozinskiǐ measure of P with respect to a matrix norm; i.e.
where I represents the identity matrix.

Theorem 4.3: The endemic equilibrium of the system (1) is globally asymptotically stable.
Proof: To show the global stability of the endemic equilibrium for the system (1), we consider a simplified case of our model by assuming φ 1 = φ 2 = α 1 = α 2 = 0; i.e. no treatments and no disease caused mortality, and apply the geometric approach summarized in Lemma 4.1. Then S + I + A + C = N is a constant which allows us to drop the last equation of system (1) and consider a three-dimensional system, written as  (23) and (24).
Note that when n = 3, the second additive compound matrix of A = (a ij ) is Thus the second additive compound matrix associated with the Jacobian is, .
According to the paper presented by Michael Y. Li [9], we need to verify the condition μ(p) < 0, which one can estimate this norm by μ(p) ≤ sup(g 1 , g 2 ). Consider that we have provide that μ < h 2 = μ + φ 2 or 0 < φ 2 which is always true. Using the equations of I * and A * , we have . This implies that μ(p) ≤ I I − μ. By the uniform persistence, there exist > 0 and T > 0 such that when t > T, we have which impliesq 2 < 0, hence, we have established the following result: From Theorem 4.3, we obtain the global asymptotically stability of the endemic equilibrium for the original system (1) under the assumptions of no treatments and disease-related mortality.

Sensitivity
There are various parameters present in our model, and the values of these parameters will directly affect the dynamics of HIV infection and transmission. These sensitivity indices are calculated by using the technique of the normalized forward sensitivity index. It will be of utmost importance to study the impact of the parameters on the basic reproduction number R 0 , determining which increase or decrease it the most, to enable implementation of the proper measures needed.

Definition 5.1:
The normalized forward sensitivity index of a variable u with respect to parameter z is given by ξ = ∂u ∂z × z u (see [3]).

Theorem 5.1:
The sensitivity and elasticity of quantity R MS with respect to the parameter β 2 , is given by K R MS β 2 = ∂R MS ∂β 2 × β 2 R MS (see [3]). Note: The transmission rates of the disease from all three infectious classes are assumed to investigate the disease dynamics.
The sensitivity indices are given in Table 1, and they are calculated by using the values of the parameters from Table 2. The positive signs of the sensitivity indices of the parameters of the basic reproduction number show that increases in the values of all of the parameters lead to an increase in the reproduction number of the disease. In Table 1, we see that the sensitivity index of β 1 , the infection rate from HIV people to Susceptibles, is 0.8856. This value signifies that R 0 is directly proportional to β 1 . Thus, the transmission rate β 1 should be controlled to reduce the number of infections and deaths of HIV/AIDS people. Other parameters give usual results such as α 1 and α 2 . These two parameters are death rates related to HIV/AIDS, so their negative signs represent normal meanings. If the values of the two parameters are increased, the reproduction number will decrease. The reduction of the reproduction number is due to the reduced number of people living with HIV/AIDS. From the sensitivity index, it is also concluded that controlling the disease transmission rates is essential to reduce the number of HIV/AIDS patients. Reducing the transmission rates will bring down the reproduction number, and it will eventually reach almost zero. As a result, there will be a small number of people living with HIV/AIDS. Thus, taking ART at early diagnosed with HIV is critical to both the health of people with HIV and the prevention of HIV transmission.

Analysis of bifurcation
We now employ the Centre Manifold theory [4] as described in Theorem 4.1 by Castillo-Chavez and Song [5], to establish the local asymptotic stability of the endemic equilibrium. Let us make the following change of variables in order to apply the Center Manifold theory: S = x 1 , I = x 2 , A = x 3 and C = x 4 , so that N = x 1 + x 2 + x 3 + x 4 . Further, by using the vector notation x = (x 1 , x 2 , x 3 , x 4 ) T , the HIV-AIDS model system (1) can be written in the The method entails evaluating the Jacobian of the system (25) at the disease-free ( 0 ) denoted J( 0 ). This gives Consider, the case when R 0 = 1. Suppose, further, that we pick β 1 = β * 1 as the bifurcation parameter. Solving for β 1 from R 0 = 1, gives It follows that J β * 1 the Jacobian of system (25) at the diseasefree with β 1 = β * 1 has a simple zero eigenvalues (with all other eigenvalues having negative real part). Hence, the Centre Manifold theory [4] can be used to analyse the dynamics of system (25).

1
For the case when R 0 = 1, it can be shown that J β * 1 has a right eigenvector (corresponding to the zero eigenvalue) given by w = (w 1 , w 2 , w 3 , w 4 ) T , where Further, the Jacobian J β * 1 has a left eigenvector (associated with the zero eigenvalue) given

Computations of a and b
It can be shown, after some algebraic manipulations (involving the associated partial nonzero partial derivative of F (at the DFE) to be used in the expression for (a) in Appendix), that Thus, a < 0 and b > 0 and Appendix item (iv) above, the following result is established: the endemic equilibrium * is locally asymptotically stable for R 0 > 1, but close to unity. A typical bifurcation diagram of I * vs. R 0 near the bifurcation point R 0 = 1 is presented in Figure 2, which is generated numerically by varying the bifurcation parameter β 1 . Figure 3 shows that all eigenvalues approach to negative real parts, thus, numerically, the endemic equilibrium point is stable.

Numerical simulation and optimal control
In this section, we will conduct a numerical simulation to verify some of our analytical results. We have run the numerical simulation for this model with various parameter values and initial conditions, and the results are consistent with our analytical predictions. In particular, using the parameter values given in Table 2 (with R 0 > 1), we display a typical set of results in Figures 3 and 4. We observe that all the curves approach the endemic equilibrium over time, indicating the global asymptotic stability of the endemic equilibrium. Next, we turn to the more general model with time-dependent controls φ 1 (t) and φ 2 (t). We consider the system on a time interval [0, T]. The functions φ 1 (t) and φ 2 (t) are assumed to be at least Lebesgue measurable on [0,T]. The control set is defined as where φ 1 max and φ 2 max denote the upper bounds for the effort of treatments. The bounds reflect practical limitation on the maximum rate of controls in a given period.
The presence of time-dependent controls makes the analysis of our system difficult. The disease dynamics now depend on the evolution of controls. In what follows we perform an optimal control study on this problem. We aim to minimize the total number of infections and the costs of control over the time interval [0, T]; i.e.
Here, the parameters c 11 , c 12 , c 21 and c 22 with appropriate units, define the appropriate costs associated with these controls. Quadratic terms are introduced to indicate nonlinear costs potentially arising at high intervention levels. The minimization process is subject to the differential equation of our system, which is now referred to as the state equations. Correspondingly, the unknowns I and A are now called the state variables, in contrast to the control variables φ 1 and φ 2 . Our goal is to determine the optimal controls φ * 1 (t) and φ * 2 (t), so as to minimize the objective functional in (27). We first establish the following theorem on the existence of optimal control. Theorem 7.1: There exist φ * 1 (t) and φ * 2 (t) ∈ such that the objective functional in (27) is minimized.
Proof: Note that the control set is closed and convex, and the integrand of the objective function in (27) is convex. Hence, based on the standard optimal control theorems in [10] the conditions for the existence of the optimal control are satisfied, as our model is linear in the control variables. Indeed, the optimal control is also unique for small T due to the Lipschitz structure of the state equations and the boundedness of the state variables.
We will follow the method described in [14] to seek the optimal control solution. This method is based on Pontryagin's Maximum Principle [25] which introduces the adjoint functions and represents an optimal control in terms of the state and adjoint functions. Essentially, this approach transfers the problem of minimizing the objective functional (under the constraint of the state equations) into minimizing the Hamiltonian concerning the controls. Our optimal control problem is as follows: Let us first define the adjoint functions λ S , λ I , λ A and λ C associated with the state equations for S, I, A and C, respectively. We then form the Hamiltonian, H, by corresponding state equations, and adding each of these products to the integrand of the objective functional. As a result, we obtain To achieve the optimal control, the adjoint functions must satisfy with transversality conditions (or final time conditions): λ S (T) = 0, λ I (T) = 0, λ A (T), and, λ C (T) = 0. The characterization of the optimal control φ * 1 (t) and φ * 2 (t) based on the condition ∂H ∂φ 1 = 0 and ∂H ∂φ 2 = 0, respectively, subject to the constraint 0 ≤ φ 1 (t) ≤ φ 1max and 0 ≤ φ 2 (t) ≤ φ 2max .
Then, 22 ). Due to the presence of both initial conditions (for the state equations) and final time conditions (for the adjoint equations), the optimal control system has to be solved numerically. We will use the Forward-Backward Sweep Method [14] to conduct the numerical simulation. First, the state equations are solved forward in time by the Euler method, with an initial guess for the control variables. Next, the adjoint equations are solved backward in time using the solutions of the state equations. The control is then updated with the new values of the state and adjoint solutions, and the process is repeated until the solutions converge.
To carry out the numerical simulation, we list the values for the various transmission rates in the state Equation (1) in Table 2. Meanwhile, we investigate the costs of ART and medical treatments for AIDS patients by setting c 11 = 1 and c 12 = 3. We further assume that for each control, the per capita cost at the quadratic level is fixed at 10, i.e.  Figures 6 and 7 show the infection curves for the model without controls (solid line), i.e. φ 1 (t) = φ 2 (t) = 0, and that with the optimal controls implemented (dashed line). We see that the number of infections for both the HIV population and AIDS patients has been reduced significantly. After the outbreak onset, the control model curve stays low due to the cooperation of ART and medical treatments. However, the results satisfy the fact that untreated HIV individuals will become AIDS patients in a period. We will observe more interesting results in our next simulation.
As can be naturally expected, the costs of these controls directly affect the strength and duration of the controls in their optimal balance. For demonstration, we now increase the cost for ART to 3 and keep the cost for medical treatment the same. Therefore, now we have a new set of cost parameters as follows: This set of cost parameters leads to some repeat cases. The reason is that with higher costs for ART, the two controls, ART and AIDS medications, are implemented with not too much average strength and shorter durations, thus some new HIV cases can occur. Based on our model formulation, we assume that AIDS patients are moved from HIV individuals  at rate σ . Thus, with a lower number of HIV infections, only a small number of individuals will move to the AIDS class. Indeed, Figures 8 and 9 show that the number of HIV individuals and AIDS patients is significantly reduced due to the incorporation of the two controls, ART and AIDS medications.
We further assume that the cost of ART is higher than that of the previous case. This time we set C 11 = 4 and keep other costs the same. Thus, now we have the set of cost parameters:   Figures 10 and 11 show the HIV and AIDS cases. We see that the number of patients is high for both. This indicates how important ART to control the disease is. Hence if one can afford to distribute ART, this control should be put into the first response to the disease outbreak.
We have also tried out various values of the transmission rate from the controlled class of individuals. The results showed that the total number of infectious people will reduce when the transmission rate from class C decreases. Therefore, reducing transmission rates from all infectious classes are important in controlling the spread of HIV.

Conclusion
In this paper, we have formulated a mathematical model of HIV-AIDS dynamics using a system of four nonlinear differential equations. We analysed the basic reproduction number, R 0 , associated with our model, and established a threshold result characterized by R 0 . If R 0 < 1, the disease will completely die out whereas if R 0 > 1, there exists a unique positive endemic equilibrium that is globally asymptotically stable. The analytical predictions are verified by our numerical simulation results for the HIV-AIDS model. Also, the focus of this study is to better understand the effects of different control strategies coupled with multiple transmission pathways of HIV, using a mathematical model that fully captures the epidemiology and pathogenesis of the disease to comprehensively explore the usefulness of public health control measures towards effective responses to HIV-AIDS epidemics in resource-constrained settings. The model consists of two types of public health controls: ART or equivalent therapy and medications for AIDS patients.
To assess the effectiveness of implementing the proposed package of public health controls, we have also performed a preliminary study based on optimal control theory to seek cost-effective solutions by assuming these interventions are time-dependent. Optimal control techniques have vastly been used to determine the effectiveness of control strategies of infectious disease dynamics. Our results suggest that different controls closely interplay with each other, and the costs of controls directly affect the length and strength of each control in an optimal strategy. We have found that due to the multiple transmission routes of HIV, a combination of multiple intervention methods achieves better results than a singlestrategy approach. Also, controlling the transmission rate from the controlled class is very important to control the spread of HIV.

Disclosure statement
No potential conflict of interest was reported by the author(s).