A TB-HIV/AIDS coinfection model and optimal control treatment

We propose a population model for TB-HIV/AIDS coinfection transmission dynamics, which considers antiretroviral therapy for HIV infection and treatments for latent and active tuberculosis. The HIV-only and TB-only sub-models are analyzed separately, as well as the TB-HIV/AIDS full model. The respective basic reproduction numbers are computed, equilibria and stability are studied. Optimal control theory is applied to the TB-HIV/AIDS model and optimal treatment strategies for co-infected individuals with HIV and TB are derived. Numerical simulations to the optimal control problem show that non intuitive measures can lead to the reduction of the number of individuals with active TB and AIDS.


1.
Introduction. According with the World Health Organization (WHO), the human immunodeficiency virus (HIV) and mycobacterium tuberculosis are the first and second cause of death from a single infectious agent, respectively [48]. Acquired immunodeficiency syndrome (AIDS) is a disease of the human immune system caused by infection with HIV. HIV is transmitted primarily via unprotected sexual intercourse, contaminated blood transfusions, hypodermic needles, and from mother to child during pregnancy, delivery, or breastfeeding [37]. There is no cure or vaccine to AIDS. However, antiretroviral (ART) treatment improves health, prolongs life, and substantially reduces the risk of HIV transmission. In both highincome and low-income countries, the life expectancy of patients infected with HIV who have access to ART is now measured in decades, and might approach that of uninfected populations in patients who receive an optimum treatment (see [12] and references cited therein). However, ART treatment still presents substantial limitations: does not fully restore health; treatment is associated with side effects; the medications are expensive; and is not curative. Following UNAIDS global report on AIDS epidemic 2013 [45], globally, an estimated 35 Mycobacterium tuberculosis is the cause of most occurrences of tuberculosis (TB) and is usually acquired via airborne infection from someone who has active TB. It typically affects the lungs (pulmonary TB) but can affect other sites as well (extrapulmonary TB). In 2012, approximately 8.6 million people fell ill with TB and 1.3 million people died from TB. Nevertheless, TB death rate dropped 45 per cent between 1990 and 2012, and 22 million lives were saved through use of strategies recommended by WHO [48].
Individuals infected with HIV are more likely to develop TB disease because of their immunodeficiency, and HIV infection is the most powerful risk factor for progression from TB infection to disease [18]. In 2012, 1.1 million of 8.6 million people who developed TB worldwide were HIV-positive. The number of people dying from HIV-associated to TB has been falling since 2003. However, there were still 320 000 deaths from HIV-associated to TB in 2012, and further efforts are needed to reduce this burden [48]. ART is a critical intervention for reducing the risk of TB morbidity and mortality among people living with HIV and, when combined with TB preventive therapy, it can have a significant impact on TB prevention [48].
Collaborative TB/HIV activities (including HIV testing, ART therapy and TB preventive measures) are crucial for the reduction of TB-HIV coinfected individuals. WHO estimates that these collaborative activities prevented 1.3 million people from dying, from 2005 to 2012. However, significant challenges remain: the reduction of tuberculosis related deaths among people living with HIV has slowed in recent years; the ART therapy is not being delivered to TB-HIV coinfected patients in the majority of the countries with the largest number of TB/HIV patients; the pace of treatment scale-up for TB/HIV patients has slowed; less than half of notified TB patients were tested for HIV in 2012; and only a small fraction of TB/HIV infected individuals received TB preventive therapy [45].
The study of the joint dynamics of TB and HIV present formidable mathematical challenges due to the fact that the models of transmission are quite distinct [36]. Some mathematical models have been proposed for TB-HIV coinfection (see, e.g., [2,3,22,28,30,36,40]). In this paper, we propose a new population model for TB-HIV/AIDS coinfection transmission dynamics, where TB, HIV and TB-HIV infected individuals have access to respective disease treatment, and single HIVinfected and TB-HIV co-infected individuals under HIV and TB/HIV treatment, respectively, stay in a chronic stage of the HIV infection.
Optimal control is a branch of mathematics developed to find optimal ways to control a dynamic system [10,16,31]. While the usefulness of optimal control theory in epidemiology is nowadays well recognized (see, e.g., [4,26,27,33,34]), and has been applied to TB models (see, e.g., [5,15,19,21,41,42]) and HIV models (see, e.g., [23,29]), to our knowledge optimal control have never been applied to a TB-HIV/AIDS coinfection model. In this paper, we apply optimal control theory to our TB-HIV/AIDS model and study optimal strategies for the minimization of the number of individuals with TB and AIDS active diseases, taking into account the costs associated to the proposed control measures.
The paper is organized as follows. The model is formulated in Section 2. The HIV-only and TB-only sub-models of the full TB-HIV/AIDS model are analyzed in Section 3 and the full TB-HIV/AIDS model is analyzed in Section 4. In Section 5 we propose an optimal control problem and apply the Pontryagin maximum principle to derive its solution. In Section 6 numerical simulations and discussion of the results are carried out for the optimal control problem associated to the TB-HIV/AIDS model. We end mentioning some possible future work in Section 7.
The susceptible population is increased by the recruitment of individuals (assumed susceptible) into the population, at a rate Λ. All individuals suffer from natural death, at a constant rate µ. Susceptible individuals acquire TB infection from individuals with active TB at a rate λ T , given by where β 1 is the effective contact rate for TB infection. Similarly, susceptible individuals acquire HIV infection, following effective contact with people infected with HIV at a rate λ H , given by where β 2 is the effective contact rate for HIV transmission. The modification parameter η A ≥ 1 accounts for the relative infectiousness of individuals with AIDS symptoms, in comparison to those infected with HIV with no AIDS symptoms. Individuals with AIDS symptoms are more infectious than HIV-infected individuals (pre-AIDS) because they have a higher viral load and there is a positive correlation between viral load and infectiousness [11]. On the other hand, η C ≤ 1 translates the partial restoration of immune function of individuals with HIV infection that use correctly ART [12].
Remark 1. For the basic and classical SIR model, one has the force of infection F given by F = βI, which models the transition rate from the compartment of susceptible individuals S to the compartment of infectious individuals I. However, for large classes of communicable diseases, it is more realistic to consider a force of infection F that does not depend on the absolute number of infectious, but on their fraction with respect to the total population N , that is, F = β I N . In our case, the force of infection for the HIV is given by λ H = β2 . Only approximately 10% of people infected with mycobacterium tuberculosis develop active TB disease. Therefore, approximately 90% of people infected remain latent. Latent infected TB people are asymptomatic and do not transmit TB [43]. Individuals leave the latent-TB class L T by becoming infectious, at a rate k 1 , or recovered, with a treatment rate τ 1 . The treatment rate for active TB-infected individuals is τ 2 . We assume that TB-recovered individuals R acquire partial immunity and the transmission rate for this class is given by β 1 λ T with β 1 ≤ 1. Individuals with active TB disease suffer induced death at a rate d T . We assume that individuals in the class R are susceptible to HIV infection at a rate λ H . On the other hand, TB-active infected individuals I T are susceptible to HIV infection, at a rate δλ H , where the modification parameter δ ≥ 1 accounts for higher probability of individuals in class I T to become HIV-positive. HIV-infected individuals (with no AIDS symptoms) progress to the AIDS class A at a rate ρ 1 , and to the class of individuals with HIV infection under treatment C H at a rate φ. Individuals in the class C H leave to the class I H at a rate ω 1 . HIV-infected individuals with AIDS symptoms are treated for HIV at the rate α 1 and suffer induced death at a rate d A . Individuals in the class I H are susceptible to TB infection at a rate ψλ T , where ψ ≥ 1 is a modification parameter traducing the fact that HIV infection is a driver of TB epidemic [24]. HIV-infected individuals (pre-AIDS) co-infected with TB-disease, in the active stage I T H , leave this class at a rate ρ 2 . A fraction p of I T H individuals take simultaneously TB and HIV treatment and a fraction q of I T H individuals take only TB treatment. Individuals in the class I T H progress to the class C H at a rate pρ 2 and to the class R H at a rate qp 2 . Individuals in the class I T H that do not take any of the TB or HIV treatments progress to the class A T at a rate (1 − (p + q))ρ 2 , and suffer TB induced death rate at a rate d T . Individuals leave L T H class at a rate τ 3 . A fraction r of L T H individuals take simultaneously TB and HIV treatment and a fraction 1 − r take only TB treatment. Individuals in the class L T H progress to the class C H at a rate rτ 3 and to the class R H at a rate (1 − r)τ 3 . Individuals in the class L T H are more likely to progress to active TB disease than individuals infected only with latent TB. In our model, this progression rate is given by k 2 . Similarly, HIV infection makes individuals more susceptible to TB reinfection when compared with non HIV-positive patients. The modification parameter associated to the TB reinfection rate, for individuals in the class R H , is given by β 2 , where β 2 ≥ 1. Individuals in this class progress to class A, at a rate ω 2 . HIV-infected individuals (with AIDS symptoms), co-infected with TB, are treated for HIV, at a rate α 2 . Individuals in the class A T suffer from AIDS-TB coinfection induced death rate, at a rate d T A . The aforementioned assumptions result in the system of differential equations Using a standard comparison theorem [25] we can show that In particular, Thus, the region Ω is positively invariant. Hence, it is sufficient to consider the dynamics of the flow generated by (3) in Ω. In this region, the model is epidemiologically and mathematically well posed [20]. Thus, every solution of the model (3) with initial conditions in Ω remains in Ω for all t > 0. This result is summarized below. 3. Analysis of the sub-models. In this section we analyze the models for HIV only (HIV-only model ) and TB only (TB-only model ).
3.1. HIV-only model. The model that considers only HIV (obtained by setting where Analogously to Lemma 2.1, we can prove that the region is positively invariant and attracting. Thus, the dynamics of the HIV-only model will be considered in Ω 1 . 3.1.1. Persistence. In this section, we look for the conditions under which the host population and disease will persist. Rewriting the submodel system (4) as in what follows we assume that β 2 (N ) is continuous for N ≥ 0 and continuously differentiable for N > 0; β 2 (N ) is monotone nondecreasing in N ; and Remark 2. In this work β 2 denotes the effective contact rate for HIV transmission. It is a constant for a concrete situation, but one can look to it as variable in the sense that, depending on the situation/region, one can have different values for this parameter. This is so because β 2 is related with the level of contagion/propagation of the disease. In Section 5 we consider fixed values for β 1 and β 2 , which represent specific cases of the infection level. By varying β 1 and β 2 we vary the basic reproduction numbers (see expressions (12) and (20) for R 1 and R 2 , respectively). Here we consider β 2 as a function of N to discuss persistence.
It is convenient to reformulate the model in terms of the fractions of the S, I H , A and C H parts of the population, and express (6) in these terms to obtain the system Equations (7) suggest that x + y + z + w = 1. The manifold x + y + z + w = 1, x, y, z, w ≥ 0, is forward invariant under the solution flow of (8), which has a global solution satisfying (7). We now show conditions under which the host population will persist.
where ε > 0 does not depend on the initial data.
Proof. We have to show that the set is uniform strong repeller for Theorem 3.2, Theorem 3.3 and Corollary 1 are taken from [3,44].
Theorem 3.2. Let X be a locally compact metric space with metric d. Let X be the disjoint union of two sets X 1 and X 2 such that X 2 is compact. Let Φ be a continuous semiflow on X 1 . Then X 2 is a uniform strong repeller for X 1 , whenever it is a uniform weak repeller for X 1 .
Theorem 3.3. Let D be a bounded interval in R and g : (t 0 , ∞) × D → R be bounded and uniformly continuous. Further, let x : (t 0 , ∞) → D be a solution of which is defined on the whole interval (t 0 , ∞). Then there exist sequences s n , t n → ∞ such that lim Corollary 1. Let the assumptions of Theorem 3.3 be satisfied. Then As the assumptions of Theorem 3.2 are satisfied, it suffices to show that X 2 is a uniform weak repeller for X 1 . Let r = y + z + w. Then, using the fact that x, y, z, w, r ≤ 1. This implies that From the equation of N in (8) we have lim inf Combining (9) and (10), we obtain that as β 2 (0) = 0 and β 2 (N ) is continuous at 0, N ∞ ≥ ε > 0 with ε not depending on the initial data. From (11) we see that we can relax β 2 (0) = 0 and require This concludes the proof.
The disease is persistent in the population if the fraction of the infected and AIDS cases is bounded away from zero. If the population dies out and the fraction of the infected and AIDS remains bounded away from zero, we would still say that the disease is persistent in the population.
Then the disease is uniformly weakly persistent insofar as with ε > 0 being independent of the initial data, provided that r(0) > 0.
The proof of Proposition 3.4 is outlined in [44].

3.1.2.
Local stability of disease-free equilibrium. The model (4) has a disease-free equilibrium (DFE), obtained by setting the right-hand sides of the equations in the model to zero, given by The linear stability of Σ 0 can be established using the next-generation operator method on the system (4). Following [46], the basic reproduction number is obtained as the spectral radius of the matrix F V −1 at the DFE, Σ 0 , with F and V given by, respectively, The basic reproduction number is given by the dominant eigenvalue of the matrix F V −1 , that is, .
The basic reproduction number R 1 represents the expected average number of new HIV infections produced by a single HIV-infected individual when in contact with a completely susceptible population [46].
Remark 3. The next-generation matrix is one of the most well known methods in epidemiology to compute the basic reproduction number for a compartmental model of the spread of infectious diseases. To calculate the basic reproduction number through this method, the whole population is divided into n compartments in which there are m < n infected compartments. Let x i , i = 1, 2, . . . , m, be the numbers of infected individuals in the ith infected compartment at time t. Now, the epidemic model is . Let x 0 denote here the disease-free equilibrium state. The Jacobian matrices of F (x) and V (x) are, respectively, where F and V are the m × m matrices given by The matrix F V −1 is known as the next-generation matrix and its spectral radius is the basic reproduction number of the model. The reader interested in all the details about the computation of the basic reproduction number by the next-generation matrix is referred to [14,46] or any good book on dynamical modeling and analysis of epidemics (e.g., [13]).
Proof. Following Theorem 2 of [46], the disease-free equilibrium, Σ 0 , is locally asymptotically stable if all the eigenvalues of the Jacobian matrix of the system (4), here denoted by M (Σ 0 ), computed at the DFE Σ 0 , have negative real parts. The Jacobian matrix of the system (4) at disease free equilibrium Σ 0 is given by One has We have just proved that the disease free equilibrium Σ 0 of model (4) is locally asymptotically stable if R 1 < 1, and unstable if R 1 > 1.
3.1.3. Global stability of disease-free equilibrium (DFE). Following [8], let us rewrite the submodel system (4) as where X = S and Z = (I H , A, C H ), with X ∈ R + denoting the total number of uninfected individuals and Z ∈ R 3 + denoting the total number of infected individuals. The disease-free equilibrium is now denoted by The conditions (H1) and (H2) below must be met to guarantee global asymptotically stability: (H1) for dX dt = F (X, 0), U 0 is globally asymptotically stable; (H2) G(X, Z) = AZ−Ĝ(X, Z),Ĝ(X, Z) ≥ 0 for (X, Z) ∈ G, where A = D Z G(U 0 , 0) is a Metzler matrix (the off diagonal elements of A are nonnegative) and G is the region where the model makes biological sense.
Setting the right hand sides of the model to zero (and noting that λ H = λ * H at equilibrium) gives Using (17) in the expression for λ * H in (16) shows that the nonzero (endemic) equilibria of the model satisfy . The force of infection at the steady-state λ * H is positive, only if R 1 > 1. We have just proved the following result.
3.1.5. Local stability of the endemic equilibrium. In what follows we prove the local asymptotic stability of the endemic equilibrium Σ H , using the center manifold theory [6], as described in [9, Theorem 4.1], with Σ H = (S * , I * H , A * , C * H ) and each of its components given as in (17). To apply this method, the following simplification and change of variables are first made. Let S = x 1 , I H = x 2 , A = x 3 and C H = x 4 , so that N = x 1 + x 2 + x 3 + x 4 . Further, by using vector notation X = (x 1 , x 2 , x 3 , x 4 ) T , the submodel (4) can be written in the form dX dt = (f 1 , f 2 , f 3 , f 4 ) T , as follows: The basic reproduction number of the submodel (4) is given by (12). Choose as bifurcation parameter β * , by solving for β 2 from R 1 = 1: The submodel (4) has a disease free equilibrium given by The Jacobian of the system (18), evaluated at Σ 0 , M (Σ 0 ), and with β 2 = β * , is given by (13). Note that the above linearized system, of the transformed system (18) with β 2 = β * , has a zero eigenvalue which is simple. Hence, the center manifold theory [6] can be used to analyze the dynamics of (18) near β 2 = β * . In particular, Theorem 4.1 in [9] is used to show the locally asymptotically stability of the endemic equilibrium point of (18), for β 2 near β * . The Jacobian M (Σ 0 ) at β 2 = β * has a right eigenvector (associated with the zero eigenvalue) given by w = [w 1 , w 2 , w 3 , w 4 ] T , where To apply Theorem 4.1 in [9] it is convenient to let f k represent the right-hand side of the kth equation of the system (18) and let x k be the state variable whose derivative is given by the kth equation for k = 1, . . . , 4. The local stability near the bifurcation point β 2 = β * is then determined by the signs of two associated constants, denoted by a and b, defined (respectively) by Note that, in f k (0, 0), the first zero corresponds to the DFE, Σ 0 , for the subsystem (4). In other words, f k (0, φ) = 0, for k = 1, . . . , 4, if and only if the right-hand sides of the equations of (4) is zero at Σ 0 . Moreover, from φ = β 2 −β * we have φ = 0 when β 2 = β * , which is the second component in f k (0, 0).
For the system (18), the associated non-zero partial derivatives at the disease free equilibrium Σ 0 are given by It follows from the above expressions that For the sign of b, it can be shown that the associated non-vanishing partial derivatives are It also follows from the above expressions that From the previous computations, we have a < 0 and b > 0. Thus, using Theorem 4.1 of [9], the following result is established.
Theorem 3.8. The endemic equilibrium Σ H is locally asymptotically stable for the basic reproduction number R 1 (12) near 1.

TB-only model. The sub-model of (3) with no HIV/AIDS disease, that is,
where The sub-model (19) was proposed and analyzed in [7]. This model incorporates the basic properties of TB transmission and dynamics. The basic reproduction number R 2 of (19) is given by The existence, uniqueness and local asymptotic stability of the disease-free and endemic equilibria are proven in [ The associated matrices F and V (see Section 3.1.2) are given, respectively, by The dominant eigenvalues of the matrix F V −1 are Thus, the basic reproduction number R 0 of the model (21) is given by Using the same procedure as in Section 3.1.2, the following result holds.

5.
Optimal control problem. In this section we present an optimal control problem, describing our goal and the restrictions of the epidemic. In the model without controls discussed so far, we have p representing the fraction of I T H individuals that take HIV and TB treatment and q representing the fraction of I T H individuals that take TB treatment only. Roughly speaking, the problem of optimal control consists to determine the optimal combination for the values of p and q. For this reason, we take p as the control u 1 and q as the control u 2 . Precisely, we add to the model (3) the two control functions u 1 (·) and u 2 (·) in the following way: As already mentioned, the controls u 1 and u 2 represent the fraction of I T H individuals that are treated for TB and HIV (simultaneously) and treated for TB only, respectively. If we consider fixed values for u 1 and u 2 in (21), then we get the model (3) with u 1 = p and u 2 = q. The aim is to find the optimal values u * 1 and u * 2 of the controls u 1 and u 2 , such that the associated state trajectories S * , L * Remark 5. Epidemiologically, our cost functional tells us that we want to minimize the number of HIV-infected individuals with AIDS symptoms co-infected with active TB. For that, one applies control measures that are associated with some implementation costs that we also intend to minimize. Other cost functionals may be used as well. Here, by considering the cost with controls in a quadratic form, we are being consistent with previous works in the literature (see, e.g., [35,42]). Moreover, a quadratic structure in the control has mathematical advantages: if the control set is a compact and convex polyhedron (as it is the case here), it imply that the Hamiltonian attains its minimum over the control set at a unique point. For future work we plan to compare the results now obtained, for a cost with a quadratic form in the controls, with those of a linear cost in the controls.
tions given in Table 1, and fixed final time T = 50 years. The initial conditions were estimated as follows. We assume that more than half of population (55%) belongs to the subgroup of susceptible and that a big percentage ( 31%) is infected with TB but is in the latent stage. This is justified from the fact that "about one-third of the world's population has latent TB", as one can find in the website of the World Health Organization (WHO) [49]. The value for the fraction of people infected with HIV is assumed 1.7%, based on HIV & AIDS Information from AVERT.org [1]: "There is either a generalised or concentrated epidemic. In a generalised epidemic, HIV prevalence is 1% or more in the general population. In a concentrated or low level epidemics, HIV prevalence is below 1% in the general population but exceeds 5% in specific at-risk populations like injecting drug users or sex workers, or HIV prevalence is not recorded at a significant level in any group." The remaining values are estimated by assuming that we are in a "controlled" situation, without large percentages in the groups of highest risk such as A, A T and C H . Our aim is to find the optimal combination of the fraction of individuals I T H that take correctly HIV and TB treatment (u * 1 ) or take only TB treatment (u * 2 ), in order to minimize the number of individuals with AIDS and TB diseases A T . Different approaches were used to obtain and confirm the numerical results. One approach consisted in using IPOPT (short for "Interior Point OPTimizer", a software library for large scale nonlinear optimization of continuous systems) [47] and the algebraic modeling language AMPL (acronym for "A Mathematical Programming Language") [17]. A second approach was to use the PROPT Matlab Optimal Control Software [32]. For more details we refer the reader to [41,42], where the same optimization approaches are used. In Figure 2 we compare the extremal dynamics C * H , I * T H and A * T associated to the extremal controls u * 1 and u * 2 with the dynamics of the model (21) with u 1 (t) = p and u 2 (t) = q, which coincide with model (3). In this simulations we consider β 1 = 0.6, β 2 = 0.1 and the rest of the parameters take the values of Table 2, which corresponds to R 0 = 4.91159 (R 1 = 4.91159, R 2 = 1.07437). We assume that the weight constants take the same value W 1 = W 2 = 50. Observe that the number of individuals with AIDS and TB diseases A T decreases significantly when the control measures u * 1 , u * 2 are implemented, see Figure 2 (c). On the other hand, the number of individuals that stays in the class C H increases in opposition to the number of individuals that have both infections HIV and TB, see Figure 2 Table 2.
of individuals I T H that take HIV and TB treatments simultaneously and only TB treatment is around 0.5 and 0.46, respectively, see Figure 3 (a) and (b). In Figure 3 (c) we observe that the extremal controls satisfy the restriction (22). At the end of 50 years, the number of TB and AIDS induced deaths reduces 5% when the controls u * 1 , u * 2 are applied, see Figure 4. Since the HIV treatments have higher costs than TB treatment, we can consider that the weight constant W 1 associated to the control u 1 takes greater values than W 2 . In this case, the fraction of individuals that take TB and HIV treatment u * 1 decreases and the fraction of individuals that take only TB treatment increases, compared to the previous case W 1 = W 2 = 50, but the associated extremal dynamics C * H , I * H and A * T behave similarly to the ones in the case W 1 = W 2 = 50, see Figure 5. The extremal controls in Figure 5 Table 2.  Figure 5. Extremal controls u * 1 and u * 2 for cost functional (23), β 1 = 0.6, β 2 = 0.1, W 1 = 500, W 2 = 50 and parameter values from Table 2. treatments is very low. If we assume that our aim is to minimize the cost functional with T = 10 years and no disease induced deaths (d T = d T A = d A = 0), that is, we wish to minimize the number of individuals that have only AIDS A and have both AIDS and TB diseases A T , the extremal controls behave in a more intuitive way. Since we assume that there is no disease induced deaths we consider that is more adequate to consider T = 10 instead of T = 50 years. Moreover, is this case, the total population is constant. We observe that the fraction of individuals that take both HIV and TB treatment u * 1 takes the maximum value for more than 7 years, and during this time the extremal control u * 2 vanishes, see Figure 6. In this case, we  Figure 6. Extremal controls u * 1 and u * 2 for cost functional J 1 , with β 1 = 0.6, β 2 = 0.1, W 1 = W 2 = 50 and parameter values from Table 2.
compare the behavior of the dynamics A, A T , I T H and C H for the following cost functionals J 1 , J 2 and J 3 , with W 1 = W 2 = 50, where that is, when both controls u 1 and u 2 are applied simultaneously, or are applied separately. The number of individuals with only AIDS, A, is lower for cost functional  J 1 and with extremal controls given by u * 1 and u * 2 in Figure 6 (a) and (b), see Figure 8 (a). However, this is not the best strategy for the reduction of the number of individuals with both AIDS and TB diseases, A T . In this case, the best strategy is to apply only control u 1 where u 1 takes the values given in Figure 7, see Figure 8 (b). The best strategy for the reduction of the total number of individuals with only AIDS and with both AIDS and TB, A + A T , during the first six years is the one associated to the cost functional J 2 , that is, apply only control u * 1 (see Figure 7 (a)) and after six years the best strategy is to apply simultaneously both controls u * 1 and u * 2 given by Figure 6. In Figure 9 we observe that the implementation of  Table 2. controls u 1 and u 2 simultaneously or separately, contribute to the reduction of the number of individuals infected with HIV and TB, I T H , and increase the number of individuals that remain in the class C H , that is, the HIV infection does not evolve to AIDS disease. In this case, the strategy of treating only TB for the individuals in the class I T H does not allow equal or better results on the reduction of individuals in the class I T H , which happens in the situation described in Figure 5 for the case W 1 = 500 and W 2 = 50.   Table 2. 7. Final comments and future work. Our numerical results only give extremals and no claims about optimality are made. As future work, it would be interesting to verify optimality by using second order optimality conditions and addressing Germany, for kindly sharing with them his expertise and for several valuable comments and helpful suggestions, which improved the quality of the paper; and to two Referees for several constructive remarks and questions.