Abstract

Tuberculosis, an airborne disease affecting almost a third of the world’s population remains one of the major public health burdens globally, and the resurgence of multidrug-resistant tuberculosis in some parts of sub-Saharan Africa calls for concern. To gain insight into its qualitative dynamics at the population level, mathematical modeling which require as inputs key demographic and epidemiological information can fill in gaps where field and lab data are not readily available. A deterministic model for the transmission dynamics of multi-drug resistant tuberculosis to assess the impact of diagnosis, treatment, and health education is formulated. The model assumes that exposed individuals develop active tuberculosis due to endogenous activation and exogenous re-infection. Treatment is offered to all infected individuals except those latently infected with multi-drug resistant tuberculosis. Qualitative analysis using the theory of dynamical systems shows that, in addition to the disease-free equilibrium, there exists a unique dominant locally asymptotically stable equilibrium corresponding to each strain. Numerical simulations suggest that, at the current level of control strategies (with Malawi as a case study), the drug-sensitive tuberculosis can be completely eliminated from the population, thereby reducing multi-drug resistant tuberculosis.

1. Introduction

Tuberculosis (TB) is a bacterial infection that is fatal if untreated timely [1]. It is an airborne disease caused by the mycobacterium tuberculosis and primarily affects the lungs (it can also affect the central nervous system, the lymphatic system, the brain, spine, and the kidneys). Approximately one-third of the world’s population is affected [2]. In 1993, concerned with the rising cases of deaths and the new infection rate which were occurring at one per second, the World Health Organization (WHO) declared TB as a global emergency. This resurgence has been closely linked with environmental and social changes that compromised people’s immune system [3]. Out of the 1.7 billion people estimated to be infected with TB, 1.3 billion lived in developing countries [2].

Active TB individuals can infect on average 10–15 other people per year if left untreated [12]. TB progression from inactive (latent) infection to active infection varies from one person to another. People suffering from AIDS have a greater risk of developing active TB with about 50% chance of developing active TB within 2 months and a 5 to 10% chance of developing active TB each year thereafter [1]. TB is treatable and curable if it is diagnosed and treated before it becomes severe [13]. WHO stresses that treatment for TB should not be undertaken unless the diagnosis is confirmed [14]. Currently five drugs are available: isoniazid, rifampicin, pyrazinamide, ethambutol, and streptomycin [13]. A combination of these drugs is required to prevent the development of drug-resistance, requiring 6–9 months of continued treatment to be effective [15].

Multidrug-resistant tuberculosis (MDR-TB) is a form of TB that is resistant to at least the two main first-line anti-TB drugs, isoniazid and rifampicin [1]. There were an estimated 0.5 million cases of MDR-TB in 2007 worldwide [14]. Drug-resistant strains are far more difficult but not impossible to treat, despite being too expensive [12]. The most important factor in preventing drug-resistant TB is to ensure full compliance with anti-TB treatment [1]. It is recommended that patients take the pills in the presence of a medical professional, an approach referred to as the directly observed therapy strategy (DOTS).

Given the scarcity of complete data, partial data obtained from the Malawi National TB Control Program [4] will be used for numerical simulations. Other parameter values are from the literature or simply assumed for the purpose of illustration. Malawi which endorsed the DOTS program since 1984 is a landlocked country in Central-Southern Africa, sharing common borders with Tanzania, Zambia, and Mozambique. The country has an estimated total population of 12.8 million and has a surface area of 118,480 km2, a quarter of which is occupied by Lake Malawi [4]. In July 2007, there was a commitment to treat all known MDR-TB cases in Malawi. By October 2007, some patients were identified, retested and a recommendation was made to start them on second-line treatment under DOTS. However, the effectiveness of the whole exercise is yet to be established as field and lab data are not yet available. Even when available, the data may not reflect the true picture because some hospitals do not collect monthly sputum specimens for checking conversion to negativity [4]. According to the 2007 tuberculosis case finding statistics, 26,299 cases were reported countrywide [4]. This is 3% less than the cases that were reported in 2006. For 2007/2008, WHO estimates that TB case detection rate for Malawi was 46%. Since TB-infected people progress faster to active TB if they are HIV positive, all TB patients are tested for HIV. Out of the 26,299 TB patients registered for anti-TB treatment, 22,744 (86%) were tested for HIV and 15,491 (68%) were found to be HIV positive [4].

Two-strain TB models that considered different interventions have been developed [5, 16, 17]. There are fundamental differences with this study. In addition to treatment, individuals are further classified based on their knowledge about health information (education) on the importance of completing their TB dosage. Also, infectious drug-sensitive individuals are diagnosed for any development of drug-resistance. Since much remains unknown about the transmission of drug-resistant TB strains, another novelty of this study is the consideration of two cases whereby an individual can get infected with MDR-TB. The first case is when latently infected individuals with drug sensitive TB come into adequate contact with an infectious MDR-TB individual and transmission takes place. The second one is when a drug sensitive TB individual can be reinfected with MDR-TB, which might also be due to incomplete treatment. Furthermore, fast and slow progression to active TB as well as endogenous re-activation and exogenous re-infection for both drug sensitive and resistant strains is accounted for.

This paper is organized as follows. In Section 2, we formulate and analyze the model. The potential impact of the various control strategies is numerically investigated in Section 3. In Section 4, we discuss the relevance of the results and possible future work.

2. Model Construction and Analysis

We consider a two-strain TB model with three interventions. The model is defined as a set of nonlinear ordinary differential equations based upon specific biological and intervention assumptions about the transmission dynamics of MDR-TB. The host population is subdivided into various classes according to their disease status: susceptible individuals (𝑆), individuals exposed to drug sensitive TB only (𝐸𝑠), infectious individuals with drug sensitive TB (𝐼𝑠), individuals exposed to MDR-TB (𝐸𝑟), infectious individuals with MDR-TB (𝐼𝑟), and individuals who have recovered from the disease (𝑅). Susceptible individuals are recruited at a constant rate, Λ. These individuals will be infected with the tubercle bacillus if they come into effective contact with an active TB case at a rate 𝜆𝑖, where the subscript 𝑖=𝑠,𝑟 denotes sensitive and MDR strains, respectively. The force of infection 𝜆𝑖 is defined as 𝜆𝑖=(𝑐𝛽𝑖𝐼𝑖)/𝑁, where 𝛽𝑖 is the probability that an individual is infected by one infectious individual, and 𝑐 is the percapita contact rate.

Progression from respective exposed classes to infectious classes is due to exogenous re-infection and endogenous reactivation. Thus, due to exogenous re-infection, individuals in 𝐸𝑠 and 𝐸𝑟 classes progress to active TB classes, 𝐼𝑠 and 𝐼𝑟, at the rate 𝛾𝑠𝜆𝑠 and 𝛾𝑟𝜆𝑟, respectively (𝛾𝑟 is the re-infection rate of exposed individuals with MDR-TB 𝛾𝑠 is similarly defined). Latently infected individuals with drug sensitive and MDR-TB strains will progress to active classes 𝐼𝑠 and 𝐼𝑟 at the rates 𝑘1 and 𝑘2, respectively, due to endogenous reactivation. Individuals in 𝐼𝑠 and 𝐼𝑟 classes are treated at the rate 𝜙𝑠 and 𝜙𝑟, respectively (realistically, it is possible that 𝜙𝑠=𝜙𝑟). They then progress to recovered class, 𝑅, if successfully treated. However, some individuals in 𝐼𝑠 class will recover naturally at a rate 𝜑 and move to 𝑅 class. Also, exposed individuals in 𝐸𝑠 and infectious individuals in 𝐼𝑠 can acquire MDR-TB if they are in contact with infectious MDR-TB individuals at a rate 𝜆𝑟 and will then enter 𝐼𝑟 class.

Infectious individuals in 𝐼𝑠 class receive treatment at a rate 𝜙𝑟, a proportion 𝑝 of which responds positively to the treatment, whereas a proportion 𝑞 partially responds to the treatment and as such they go back to 𝐸𝑠 class. The remaining proportion (1(𝑝+𝑞)) will not complete the treatment which may result in the development of MDR-TB and these individuals move to 𝐸𝑟 class. In addition, health education is offered to infectious individuals with drug sensitive strains only at a rate 𝑎. This is due to the nature of the disease, that is, one is diagnosed with drug sensitive TB (at a rate 𝜎 in this case) which later progress to MDR-TB if treatment compliance is disregarded [13]. Both 𝜙𝑟 and 𝜎 also describe a consequence of incomplete treatment, and as such, treatment rate 𝜙𝑟 is also a result of a diagnosis.

Susceptible individuals who become infected progress faster to active drug sensitive TB, that is, from 𝑆 to 𝐼𝑠 class at a rate 𝜌𝑠 and to resistant strain class 𝐼𝑟, at a rate 𝜌𝑟; this might be due to other immunocompromised factors such as HIV and malaria that weakens individuals’ immune systems leaving them very vulnerable to TB attack. Thus, (1𝜌𝑠) and (1𝜌𝑟) denote slow progression to active drug sensitive and MDR strains, respectively. We assume that recovery is non permanent and as such recovered individuals are infected with drug sensitive TB at a rate 𝜆𝑠, move to 𝐸𝑠 class where they become infected with MDR-TB at a rate 𝜆𝑟 to move into the 𝐸𝑟 class. Furthermore, infectious individuals in 𝐼𝑠 class die due to the disease at a rate 𝑑𝑠 and those in 𝐼𝑟 class die at a rate 𝑑𝑟. All individuals in different subgroups die naturally at a rate 𝜇. A schematic diagram of the model is depicted in Figure 1, and the associated parameters are described in Table 1.

With the pervious assumptions, terminology and inter-relations between the parameters and variables as described by Figure 1, the dynamics of the MDR-TB model can be described by the following deterministic system of nonlinear ordinary differential equations: 𝑆𝜆(𝑡)=Λ𝑠+𝜆𝑟𝐸𝑆𝜇𝑆,𝑠(𝑡)=1𝜌𝑠𝜆𝑆+𝑅𝑠𝛾𝑠𝜆𝑠+𝜆𝑟𝐸𝑠𝜙𝑠+𝑘1𝐸+𝜇𝑠+𝑞𝜙𝑟𝐼𝑠,𝐼𝑠(𝑡)=𝜌𝑠𝜆𝑠𝛾𝑆+𝑠𝜆𝑠+𝑘1𝐸𝑠𝜆𝑟𝐼𝑠𝑑𝑠+𝑎𝑝𝜙𝑟+𝜑𝑠𝐼+𝜎+𝜇𝑠,𝐸𝑟(𝑡)=1𝜌𝑟𝜆𝑆+𝑅𝑟+(1(𝑝+𝑞))𝜙𝑟𝐼𝑠𝜆𝑠+𝛾𝑟𝜆𝑟+𝑘2𝐸+𝜇𝑟,𝐼𝑟(𝑡)=𝜌𝑟𝜆𝑟𝜆𝑆+𝑠+𝛾𝑟𝜆𝑟+𝑘2𝐸𝑟+𝜆𝑟𝐸𝑠+𝜆𝑟𝐼𝑠+𝜎𝐼𝑠𝑑𝑟+𝜙+𝜑𝑟𝐼+𝜇𝑟,𝑅(𝑡)=𝑎𝑝𝜙𝑟+𝜑𝑠𝐼𝑠+𝜙𝑠𝐸𝑠+𝜙+𝜑𝑟𝐼𝑟𝜆𝑠+𝜆𝑟𝑅𝜇𝑅,(1) where the force of infection 𝜆𝑠=𝑐𝛽𝑠(𝐼𝑠/𝑁), 𝜆𝑟=𝑐𝛽𝑟(𝐼𝑟/𝑁). The initial conditions are 𝑆(0)=𝑆0,𝐸𝑠(0)=𝐸0𝑠,𝐼𝑠(0)=𝐼0𝑠,𝐸𝑟(0)=𝐸0𝑟,𝐼𝑟(0)=𝐼0𝑟,𝑅(0)=𝑅0. The total population 𝑁 (say) of system (1) is given by 𝑁=𝑆+𝐸𝑠+𝐼𝑠+𝐸𝑟+𝐼𝑟+𝑅. Model system (1) monitors a human population; therefore, all its associated parameters and state variables are assumed to be nonnegative for all 𝑡>0. Thus, the feasible solutions of system (1) are well-defined in Γ=𝑆(𝑡),𝐸𝑠(𝑡),𝐼𝑠(𝑡),𝐸𝑟(𝑡),𝐼𝑟(𝑡),𝑅(𝑡)6+Λ𝑁𝜇,(2) which is positively invariant and attracting and it is sufficient to consider solutions in Γ [18]. Furthermore, existence, uniqueness, and continuation of results for system (1) hold in this region. Also, all solutions of model system (1) starting in Γ remain in Γ for all 𝑡0.

2.1. The Disease-Free Equilibrium and Its Stability

In the absence of infection (i.e., 𝐸𝑠=𝐸𝑟=𝐼𝑠=𝐼𝑟=0), model system (1) has a disease-free equilibrium 𝐸0 given by 𝐸0=𝑆0,𝐸0𝑠,𝐼0𝑠,𝐸0𝑟,𝐼0𝑟,𝑅0=Λ𝜇,0,0,0,0,0.(3) The potential intensity of transmission and the dynamics of a disease are often investigated in terms of the reproductive number, which represents the mean number of secondary cases a typical single infected individual will generate in a totally naive/susceptible population during his/her entire period of infectiousness. The linear stability of the disease-free equilibrium 𝐸0 is investigated using the next generation matrix for system (1) [19]. To this effect, we compute the effective reproduction number 𝑅𝑒, the threshold for endemic persistence and epidemic spread of the disease. This is an important nondimensional quantity in epidemiology as it sets the threshold for predicting a disease outbreak and for evaluating its control strategies [20]. Therefore, whether a disease becomes persistent or dies out in a community depends on the size of this threshold parameter. Mathematically, 𝑅𝑒 is the spectral radius of the next-generation matrix [19]. The next-generation matrix calculation (see details in Appendix A) shows that the effective reproduction number (or epidemic threshold) is 𝑅𝑒𝑅=max𝑠,𝑅𝑟,(4) where 𝑅𝑠=𝑐𝛽𝑠𝜇𝜌𝑠+𝜙𝑠𝜌𝑠+𝑘1𝜙𝑠+𝑘1𝑑+𝜇𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇𝑞𝜙𝑟𝑘1,𝑅𝑟=𝑐𝛽𝑟𝑘2+𝜇𝜌𝑟𝑘2𝑑+𝜇𝑟+𝜙+𝜑𝑟.+𝜇(5)𝑅𝑠 and 𝑅𝑟 are, respectively, the reproduction numbers for drug-sensitive TB strain only and MDR-TB strain only. 𝑅𝑒 measures the average number of new infections generated by a typical infectious individual in a community where intervention strategies are in place. Thus, in the absence of diagnosis, treatment, and health education (i.e., 𝜙𝑠=𝜙𝑟=𝜙=𝑎=𝜎=0), (A.7) reduces to 𝑅0𝑠=𝑐𝛽𝑠𝑘1+𝜇𝜌𝑠𝑘1𝑑+𝜇𝑠+𝜑𝑠,𝑅+𝜎+𝜇0𝑟=𝑐𝛽𝑟𝑘2+𝜇𝜌𝑟𝑘2𝑑+𝜇𝑟+𝜑𝑟,+𝜎+𝜇(6)𝑅0=max{𝑅0𝑠,𝑅0𝑟}. The threshold quantity 𝑅0 is the basic reproduction number of infection representing the average number of new infections generated by a single infective individual in a completely naive population. Each term in 𝑅𝑠 and 𝑅𝑟 has an epidemiological interpretation. For the drug-sensitive reproduction number,(i)𝑘1/(𝜙𝑠+𝑘1+𝜇) is the expected fraction of individuals that will progress from 𝐸𝑠 class to 𝐼𝑠;(ii)1/(𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇) is the expected time infectious individuals with drug-sensitive TB spend in 𝐼𝑠 class.

A similar interpretation caters for the drug-resistant reproduction number. Thus, from [19] the following result holds.

Theorem 1. The disease-free equilibrium 𝐸0 of model system (1) is locally asymptotically stable if 𝑅𝑒<1, that is, 𝑅𝑠<1 and 𝑅𝑟<1, and unstable if 𝑅𝑒>1, that is, 𝑅𝑠>1 and 𝑅𝑟>1.

2.2. The Endemic Equilibria

For system (1), there are three possible endemic equilibria; two boundary equilibrium points which are 𝐸1 (exists only when drug-sensitive strain is present) and 𝐸2 (exists only when drug-resistant strain is present) and the equilibrium point 𝐸3 which exists when both strains are present or coexist.

2.2.1. The Drug-Sensitive TB-Only Endemic Equilibrium

This is obtained by setting classes 𝐸𝑟=𝐼𝑟=0. This reduces system (1) to 𝑆(𝑡)=Λ𝜆𝑠𝐸𝑆𝜇𝑆,𝑠(𝑡)=1𝜌𝑠𝜆𝑆+𝑅𝑠𝛾𝑠𝜆𝑠𝐸𝑠𝜙𝑠+𝑘1𝐸+𝜇𝑠+𝑞𝜙𝑟𝐼𝑠,𝐼𝑠(𝑡)=𝜌𝑠𝜆𝑠𝛾𝑆+𝑠𝜆𝑠+𝑘1𝐸𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠𝐼+𝜎+𝜇𝑠,𝑅(𝑡)=𝑎𝑝𝜙𝑟+𝜑𝑠𝐼𝑠+𝜙𝑠𝐸𝑠𝜆𝑠+𝜇𝑅.(7) The drug-sensitive TB-only equilibrium in terms of the equilibrium value of the force of infection 𝜆𝑠 is given by 𝐸1=(𝑆,𝐸𝑠,𝐼𝑠,𝑅,0,0), where 𝑆=Λ𝜆𝑠+𝜇,𝐸𝑠=𝑎1𝜆𝑠2+𝑎2𝜆𝑠𝜆𝑠𝑏+𝜇1𝜆𝑠2+𝑏2𝜆𝑠+𝑏3,𝐼𝑠=𝜌𝑠𝜆𝑠Λ𝑏1𝜆𝑠2+𝑏2𝜆𝑠+𝑏3+𝛾𝑠𝜆𝑠+𝑘1𝑎1𝜆𝑠2+𝑎2𝜆𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠𝜆+𝜎+𝜇𝑠𝑏+𝜇1𝜆𝑠2+𝑏2𝜆𝑠+𝑏3,𝑅=𝜆𝑠𝑎3𝜆𝑠2+𝑎4𝜆𝑠+𝑎5𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠𝜆+𝜎+𝜇𝑠𝑏+𝜇1𝜆𝑠2+𝑏2𝜆𝑠+𝑏3,(8) with 𝑎1=Λ1𝜌𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑞𝜙𝑟𝜌𝑠𝑎𝑝𝜙𝑟+𝜑𝑠𝜌𝑠,𝑎2=𝜇Λ1𝜌𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑞𝜙𝑟𝜌𝑠,𝑎3=𝑏1𝑎𝑝𝜙𝑟+𝜑𝑠Λ𝜌𝑠+𝑎1𝑎𝑝𝜙𝑟𝛾𝑠,𝑎4=𝑏2𝑎𝑝𝜙𝑟+𝜑𝑠Λ𝜌𝑠+𝑎1𝑎𝑝𝜙𝑟+𝜑𝑠𝑘1+𝜙𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑎2𝑎𝑝𝜙𝑟+𝜑𝑠𝛾𝑠,𝑎5=𝑏3𝑎𝑝𝜙𝑟+𝜑𝑠Λ𝜌𝑠+𝑎2𝑎𝑝𝜙𝑟+𝜑𝑠𝑘1+𝜙𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠,𝑏+𝜎+𝜇1=𝛾𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇𝛾𝑠𝑞𝜙𝑟𝑎𝑝𝜙𝑟+𝜑𝑠𝛾𝑠,𝑏2=𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠𝜙+𝜎+𝜇𝑠+𝑘1+𝜇+𝜇𝛾𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇𝑞𝜙𝑟𝜇𝛾𝑠+𝑘1𝑎𝑝𝜙𝑟+𝜑𝑠𝑘1+𝜙𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠,𝑏+𝜎+𝜇3𝑑=𝜇𝑠+𝑎𝜙𝑟+𝜑𝑠𝜙+𝜎+𝜇𝑠+𝑘1+𝜇𝜇𝑞𝜙𝑟𝑘1.(9) Substituting 𝐸1 into the relationship 𝜆𝑠=(𝑐𝛽𝑠𝐼𝑠)/𝑁, we obtain the drug-sensitive TB-only endemic equilibrium that satisfies the following polynomial: 𝜆𝑠𝜆𝑠=𝜆𝑠𝐴1𝜆𝑠2+𝐵1𝜆𝑠+𝐶1=0,(10) where 𝐴1=𝑎1𝜇𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑘1+𝑎2𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝛾𝑠+𝑏2Λ𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝜌𝑠+𝑏1𝑑𝜇Λ𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑎4𝑐𝛽𝑠𝜌𝑠Λ𝑏1+𝑎1𝛾𝑠,𝐵1=𝑎2𝑘1𝑑+𝜇𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑏3Λ𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝜌𝑠+𝑏2𝑑Λ𝜇𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝑎5𝑐𝛽𝑠𝑏3𝜌𝑠Λ+𝑎2𝑘1,𝐶1=𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇+𝜌𝑠𝜙𝑠+𝑘1+𝜇𝑘1𝑞𝜙𝑟×1𝑅𝑠.(11) The solution 𝜆𝑠=0 in (10) corresponds to the disease-free equilibrium and (𝜆𝑠)=0 corresponds to the existence of multiple equilibria. The coefficient 𝐴1 is always positive and 𝐶1 is positive if 𝑅𝑠 is less than unity and negative if 𝑅𝑠 is greater than unity. Thus, we have established the following result.

Theorem 2. The drug sensitive TB only model system (7) has(i)precisely one unique endemic equilibrium if 𝐶1<0𝑅𝑠>1,(ii)precisely one unique endemic equilibrium if 𝐵1<0 and 𝐶1=0 or 𝐵214𝐴1𝐶1=0,(iii)precisely two endemic equilibria if 𝐶1>0, 𝐵1<0 and 𝐵214𝐴1𝐶1>0,(iv)otherwise, there is no endemic equilibrium.

Condition (iii) implies that the dynamical phenomenon of backward bifurcation where a stable endemic equilibrium coexists with a stable disease-free equilibrium when the associated reproduction number is less than unity. This has important implications for disease control. In such a scenario, the classical requirement of the reproduction number being less than unity becomes only a necessary, but not a sufficient condition for disease elimination. To find the backward bifurcation point, we set the discriminant 𝐵214𝐴1𝐶1 to zero. Making 𝑅𝑠 the subject of formula, we obtain 𝑅𝑐𝑠𝐵=1214𝐴1𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠𝜙+𝜎+𝜇𝑠+𝑘1+𝜇𝑘1𝑞𝜙𝑟.(12) Hence, it can be shown that backward bifurcation occurs for values of 𝑅𝑠 in the range 𝑅𝑐𝑠<𝑅𝑠<1 (see Figure 2). Figure 2 shows a backward bifurcation diagram of model system (7). From the diagram, we observe that if 𝑅𝑠<1 and then increases to unity, the number of TB cases increases abruptly thus, the disease-free equilibrium coexist with the endemic equilibrium implying that; the disease can invade the population to a relatively high endemic level. In addition, decreasing 𝑅𝑠 to its former level will not necessarily make the disease disappear. This is a consequence of the exogenous re-infection feature of TB. Hence, exogenous re-infection is capable of sustaining TB even when the reproduction number is below one [21]. However, backward bifurcation is illustrated by specific choice of parameters which may not be epidemiologically realistic.

2.2.2. The Drug-Resistant TB Only Endemic Equilibrium

This is obtained by setting 𝐸𝑠=𝐼𝑠=0 in model system (1). Hence, system (1) becomes 𝑆(𝑡)=Λ𝜆𝑟𝐸𝑆𝜇𝑆,𝑟(𝑡)=1𝜌𝑟𝜆𝑟𝑆+𝜆𝑟𝛾𝑅𝑟𝜆𝑟+𝑘2𝐸+𝜇𝑟,𝐼𝑟(𝑡)=𝜌𝑟𝜆𝑟𝛾𝑆+𝑟𝜆𝑟+𝑘2𝐸𝑟𝑑𝑟+𝜙+𝜑𝑟𝐼+𝜇𝑟,𝑅(𝑡)=𝜙+𝜑𝑟𝐼𝑟𝜆𝑟+𝜇𝑅,(13) so that 𝑁=𝑆+𝐸𝑟+𝐼𝑟+𝑅. Therefore, the drug-resistant TB only equilibrium in terms of the equilibrium value of the force of infection 𝜆𝑟 is given by 𝐸2=(𝑆,0,0,𝐸𝑟,𝐼𝑟,𝑅)=(𝑆,𝐸𝑟,𝐼𝑟,𝑅), where 𝑆=Λ𝜆𝑟,𝐸+𝜇𝑟=1𝜌𝑟Λ𝜆𝑟𝑚4𝜆𝑟2+𝑚5𝜆𝑟+𝑚6𝛾𝑟𝜆𝑟+𝑘2𝜆+𝜇𝑟𝑚+𝜇4𝜆𝑟2+𝑚5𝜆𝑟+𝑚6+𝜙𝜆𝑟𝑚1𝜆𝑟2+𝑚2𝜆𝑟+𝑚3𝛾𝑟𝜆𝑟+𝑘2𝜆+𝜇𝑟𝑚+𝜇4𝜆𝑟2+𝑚5𝜆𝑟+𝑚6,𝐼𝑟=𝑚1𝜆𝑟2+𝑚2𝜆𝑟+𝑚3𝑚4𝜆𝑟2+𝑚5𝜆𝑟+𝑚6,𝑅=𝜙𝑚1𝜆𝑟2+𝑚2𝜆𝑟+𝑚3𝜆𝑟𝑚+𝜇4𝜆𝑟2+𝑚5𝜆𝑟+𝑚6,(14) with 𝑚1=Λ𝛾𝑟,𝑚2=𝜇Λ𝜌𝑟,𝑚3=Λ𝑘2,𝑚4=𝛾𝑟𝑑𝑟,𝑚+𝜇5=𝑑𝑟+𝜙+𝜑𝑟𝑘+𝜇2+𝜇+𝜇𝛾𝑟𝜙+𝜑𝑟𝑘2,𝑚6𝑘=𝜇2𝑑+𝜇𝑟+𝜙+𝜑𝑟.+𝜇(15) Substituting 𝐸2 into the equation 𝜆𝑟=(𝑐𝛽𝑟𝐼𝑟)/𝑁, we obtain an endemic equilibrium of the drug-resistant TB only that satisfies the polynomial given by 𝜆𝑟𝑔𝜆𝑟=𝜆𝑟𝐴𝜆𝑟2+𝐵𝜆𝑟+𝐶=0,(16) where 𝐴=𝑚4Λ𝑘2+𝜇+𝑚5Λ𝛾𝑟+1𝜌𝑟+𝑚1𝑘2+𝜇𝜙+𝜑𝑟+𝜇+𝑚2𝜙+𝜑𝑟+𝛾𝑟𝜇+𝑘2+𝜇+𝜙+𝜑𝑟𝛾𝑟+𝑚3𝛾𝑟𝑐𝛽𝑟𝜇𝛾𝑟+𝑘2𝑚+𝜇1+𝑚2𝛾𝑟,𝐵=𝑚5Λ𝑘2+𝜇+𝑚6Λ𝛾𝑟+1𝜌𝑟+𝑚2𝑘2+𝜇𝜙+𝜑𝑟+𝜇+𝑚3𝜙+𝜇𝛾𝑟+𝑘2+𝜇+𝜙+𝜑𝑟𝛾𝑟𝑐𝛽𝑟𝑚1𝑘2+𝜇𝜇+𝑚2𝜇𝛾𝑟+𝑘2+𝜇+𝑚3𝛾𝑟,𝑘𝐶=𝜇2𝑑+𝜇𝑟+𝜙+𝜑𝑟+𝜇𝜇𝛾𝑟+𝑘2+𝜇1𝑅𝑟.(17) The root 𝜆𝑟=0 in (16) corresponds to 𝐸0 (its stability has already been established) and 𝑔(𝜆𝑟)=0 can be analyzed for the possibility of existence of multiple equilibria. It is worth mentioning here that the coefficient 𝐴 is always positive and 𝐶 is positive if 𝑅𝑟<1 and negative if 𝑅𝑟>1, hence, the following result.

Theorem 3. The drug-resistant TB only model (13) has(1)precisely one unique endemic equilibrium if 𝐶<0𝑅𝑟>1,(2)precisely one unique endemic equilibrium if 𝐵<0 and 𝐶=0 or 𝐵24𝐴𝐶=0,(3)precisely two endemic equilibria if 𝐶>0, 𝐵<0 and 𝐵24𝐴𝐶>0,(4)otherwise there is no endemic equilibrium.

The backward bifurcation point can be found by setting the discriminant 𝐵24𝐴𝐶 to zero. Then, making 𝑅𝑟 the subject of the formula, we obtain 𝑅𝑐𝑟𝐵=12𝑘4𝐴𝜇2𝑑+𝜇𝑟+𝜙+𝜑𝑟+𝜇𝜇𝛾𝑟+𝑘2,+𝜇(18) from which it can be shown that backward bifurcation occurs for values of 𝑅𝑟 in the range 𝑅𝑐𝑟<𝑅𝑟<1. The following result provides a condition for the existence of the drug-resistant TB only endemic equilibrium point, 𝐸2.

Theorem 4. The drug-resistant TB only endemic equilibrium 𝐸2 exists whenever 𝑅𝑟>1.

Proof. Solving for 𝜆𝑟 in (16), we obtain 𝜆𝑟=𝐵+𝐵24𝐴𝐶.2𝐴(19) The disease is endemic when 𝜆𝑟>0 which occurs if and only if 𝐵24𝐴𝐶>𝐵2𝑘4𝐴𝜇2𝑑+𝜇𝑟+𝜙+𝜑𝑟×+𝜇𝜇𝛾𝑟+𝑘2+𝜇1𝑅𝑟<0𝑅𝑟>1.(20)
Thus, 𝐸2 exists whenever 𝑅𝑟>1.

Again, using the Center Manifold theory [22], the local asymptotic stability of 𝐸2 is established (see details in Appendix B). The bifurcation diagram of the drug-resistant TB only model is illustrated in Figure 3.

Figure 3 illustrates a case of a backward bifurcation of system (13). As 𝑅𝑟 approaches unity, it can be seen from the diagram that the number of secondary transmission suddenly increases giving rise to a situation whereby the disease-free equilibrium coexist with the endemic equilibrium.

2.2.3. Two-Strain Model: Drug Sensitive and MDR-TB Endemic Equilibrium

Having analyzed the dynamics of the two submodels, the full drug sensitive and MDR-TB model is now considered. Its endemic equilibrium occurs when both drug sensitive and MDR-TB strains circulate in the community and is denoted by 𝐸3=𝑆,𝐸𝑠,𝐼𝑠,𝐸𝑟,𝐼𝑟,𝑅.(21) It is a daunting task to explicitly express 𝐸3 in terms of the equilibrium value of the force of infection 𝜆𝑟𝑠. As in the previous sections, it can be shown, using the next-generation method that the associated reproduction number for the full model is 𝑅𝑒=max{𝑅𝑠,𝑅𝑟}, where 𝑅𝑠 and 𝑅𝑟 are, respectively, the reproduction numbers of drug sensitive and drug-resistant TB only sub-models given earlier. 𝑅𝑒=max{𝑅𝑠,𝑅𝑟} implies that the two TB strains (drug sensitive and drug-resistant) escalate each other and competitive exclusion may occur.

If 𝑅𝑒=max{𝑅𝑠,𝑅𝑟}, then from Theorem 2, the drug sensitive TB only sub-model has a backward bifurcation for values of 𝑅𝑠 such that 𝑅𝑐𝑠<𝑅𝑠<1 and Theorem 3 showed that the drug-resistant TB only sub-model exhibits backward bifurcation for values of 𝑅𝑟 such that 𝑅𝑐𝑟<𝑅𝑟<1. Thus, the two-strain model will also exhibit the phenomenon of backward bifurcation whenever 𝑅𝑒<1, and consequently, the coexistence endemic equilibrium 𝐸3 is only locally asymptotically stable when 𝑅𝑒>1.

The existence of multiple endemic equilibria is of general interest far beyond tuberculosis epidemiology. An important principle in theoretical biology is that of competitive exclusion which states that no two species can forever occupy the same ecological niche [23]. The system studied has the requisite in which the principle of competitive exclusion holds. Since model (1) exhibits the phenomenon of backward bifurcation thereby implying that the two-strain model is only locally stable, the strain with the large reproduction number colonizes the other strain [24].

3. Model Simulations

Graphical representations to support the analytical results are provided using a set of parameter values given in Table 1. These values were obtained from the National Tuberculosis Control Programme secretariat (Lilongwe, Malawi). Incomplete data from the National TB Control Program proves to be a major challenge, and the actual values of most of the parameters are not readily available [25]. Therefore, we use values from the literature for the purpose of illustration. We simulate both the drug sensitive and MDR-TB dynamics in the absence of any intervention and when the interventions are present as well as the effect of varying each intervention parameter on the number of infected populations. All figures are generated using the parameter values presented in Table 1 and the following initial conditions 𝑆0=14000,𝐸0𝑠=10500,𝐼0𝑠=7500,𝐸0𝑟=6500,𝐼0𝑟=5500,𝑅=4000. The rationale behind this particular choice of the initial conditions is to capture the dynamics of the epidemic in a small community. TB is a disease with slow dynamics and consequently, TB epidemics must be studied and assessed over extremely long windows in time [26]. The time unit throughout is per year.

(a) Absence of any Intervention Strategy
In the absence of interventions, the susceptible population initially decreases and then increases to its carrying capacity as shown by Figure 4(a). On the other hand, the populations of latent drug sensitive TB, infectious drug sensitive TB, latent drug-resistant TB, and infectious drug-resistant TB decrease to an endemic level with increasing time as shown by Figures 4(b), 4(c), 4(d), and 4(e), respectively. This indicates that as long as there are no interventions to control the spread of the disease, the disease will not clear from the population since the basic reproduction number 𝑅0=1.4286>1. This result supports the theorem on local stability of endemic equilibrium.

(b) With Control Strategies (Presence of Interventions)
When interventions are introduced, improved trends of the populations are observed. For instance, in Figure 5(a), the susceptible population increases compared to the case when no interventions are available. Further improved trends can also be seen in Figures 5(b)5(f). Figures 5(b) and 5(c) indicate that individuals infected with drug sensitive TB decrease and eventually diminish to zero as a result of the interventions. This means that, if the disease threshold is below unity (𝑅𝑒=0.2987), both drug sensitive and resistant strains can be eliminated. Figure 5(e) depicts the time series plot of the population density of infectious individuals with drug-resistant TB which decreases but does not tend to zero. This simply means that the interventions in place are not enough to completely eradicate the epidemic from the population. The observed decrease of the number of drug-resistant TB individuals may be the result of abrupt reductions in the rates of disease progression [27].

MDR-TB which is difficult to treat, spreads fastest in areas with poor adherence to second-line drug. This poor adherence is frequently caused by shortages of second-line drugs which are quite expensive and as such minimal treatment is offered to those infected. Figure 5(f) shows that the recovered population increases as a result of the interventions unlike in Figure 4(f). In other words, we observe that the introduction of treatment, diagnosis, and health education in a TB stricken community reduces the impact of MDR-TB but cannot completely clear it from the community, because higher levels of treatment may lead to an increase of the epidemic size, and the extend to which this occurs depends on the factors such as drug efficacy and resistance development [28]. Figure 6(a) shows that diagnosis of infectious individuals with drug sensitive TB to determine whether or not the infection has developed resistance to drugs is very crucial in MDR-TB control. More infectious individuals with sensitive TB needs effective treatment that should correlate with the diagnosis levels. In addition, diagnosis is very important to detect the number of people who have developed resistant strains and are eligible to the second-line treatment to prevent the infection from a possible spread. As for the sick individuals with MDR-TB, treatment of the infection is paramount as indicated by Figure 6(b). Also, from Figure 6(b), education campaigns alone cannot curb or reduce the infection but work hand in hand with treatment as well as diagnosis. In other words, Figure 6 suggests that all individuals diagnosed with MDR-TB should be educated on the importance of treatment compliance and completion.

3.1. Dynamics of the Populations under Different Interventions
3.1.1. The Effect of Treatment on MDR-TB Dynamics

It is assumed under this scenario that treatment is given to latent and infectious individuals with drug sensitive TB as well as infectious individuals with MDR-TB. We then investigate the impact of each of these control measures on all the infected populations of both strains. As the treatment rate of latent individuals with drug sensitive TB, 𝜙𝑠, increases, 𝑅𝑒 decreases so are the infectious populations with both strains. Thus, treating more latent individuals with drug sensitive TB can eliminate drug sensitive TB (Figures 10(a) and 10(c)). This is the case because as more latent individuals with drug sensitive TB are treated, then only a few of them will progress to active infection. Also, increasing 𝜙𝑠 reduces the number of infectious individuals with MDR-TB since the treatment will prevent the infection from developing resistance to the drugs. Although, this is the case, MDR-TB may not completely be eradicated from the population due to re-infection and relapse as illustrated in Figure 10(d), and also due to the fact that treatment efficacy is less than 100%. Figures 7(a) and 7(b) show that increasing 𝜙𝑟 reduces both 𝑅𝑒 and the latent and infectious populations with drug sensitive TB to zero over time. Treating more infectious individuals with drug sensitive TB prevents the infection from developing resistance to drugs and hence reduces the number of infectious population with MDR-TB as shown in Figure 7(d). However, Figure 7(d) also shows that, at this level of treatment, MDR-TB cannot be absolutely wiped out of the society which confirms the complexity of the disease. Figures 8(a) and 8(b) show that as the treatment rate of infectious individuals with MDR-TB, 𝜙, increases, 𝑅𝑟 reduces to less than unity and decreasing trends for latent and infectious individuals with MDR-TB are observed although they do not decay to zero due to the continuous development of resistance as treatment is not fully (or 100%) effective.

3.1.2. The Effect of Diagnosis on MDR-TB Dynamics

Figure 9 shows that increasing the value of 𝜎 reduces 𝑅𝑒 and also decreases the infectious populations with drug sensitive and MDR-TB. From Figure 9(a), drug sensitive TB can be completely eliminated from the population if more people are diagnosed. This is mainly the case because usually diagnosis leads to treatment which reduces the infection (Figures 10, 7, and 8). On the contrary, diagnosis of more infectious individuals with drug sensitive TB is not a guarantee of eliminating MDR-TB as it only reduces the number of infectious individuals with MDR-TB but does not wipe the infection out of the population as illustrated by Figure 9(b). Therefore, increase in diagnosis should be correlated with increase in treatment to ensure treatment for all infectious individuals after they are detected.

3.1.3. The Effect of Health Education on MDR-TB Dynamics

Figure 11 illustrates the importance of health education in the fight against MDR-TB. It is observed in Figure 11(a) that when more people receive health education on the importance of adhering to the doctor’s recommendation on how to take their TB regimens, the infectious population with drug sensitive TB decreases and eventually decays to zero. Also, this strategy reduces 𝑅𝑒 to further smaller values. Consequently, health education slightly reduces infectious individuals with MDR-TB as shown in Figure 11(b). This is possible because treatment adherence and compliance reduce the likelihood of the infection developing drug-resistance. However, Figure 11(b) also indicates that education alone is not enough to completely eliminate MDR-TB from the community as not all people will follow these rational instructions. In addition, exogenous re-infection and regression also make the efforts to root out MDR-TB difficult but not impossible. Thus, preventing re-infection and regression are viable. Figure 12 shows that as more infectious individuals with drug sensitive TB receive health education, the number of recovered population increases. This supports the fact that education has a positive impact on TB dynamics as depicted in Figure 11. Thus, educating more infectious individuals with drug sensitive TB increases the number of people recovering from the infection which is a positive development for the management of MDR-TB. Therefore, stepping up TB information/awareness campaigns should be given prominence in TB control programmes.

4. Discussion and Conclusion

A two-strain TB model with diagnosis, treatment, and health education is formulated and analyzed. The main objective of this theoretical study was to assess the impact of these control strategies on the transmission dynamics of MDR-TB (with Malawi as a case study). We note however that the results presented are general and can be applied to other settings because neither the model, nor the parameters values represent characteristics unique of Malawi. The effective reproduction number was computed and used to compare the effect of each intervention strategy on the MDR-TB dynamics.

Using the theory of dynamical systems, qualitative analysis shows that the model has two equilibria the disease-free equilibrium and endemic equilibrium. Using the next-generation operator approach, it was found that, whenever the threshold that describes endemic persistence of the disease, 𝑅𝑒<1 (i.e., 𝑅𝑠<1 and 𝑅𝑟<1), the disease-free equilibrium is locally asymptotically stable and becomes unstable whenever 𝑅𝑒>1 (𝑅𝑠>1 and 𝑅𝑟>1). The existence and stability of the endemic equilibrium was determined using the Center Manifold theory [22]. Near the threshold 𝑅𝑒=1, there exists a stable endemic equilibrium which is locally asymptotically stable for 𝑅𝑒>1. In the absence of interventions, the effective reproduction number, 𝑅𝑒, reduces to the basic reproduction number 𝑅0. As customary in epidemiological models, the disease-free and endemic equilibria are found and their stability is investigated depending on the system parameters. Because of the occurrence of backward bifurcation in some parameter regimes, the system exhibits a bistability between a disease-free and endemic steady states. Whether the parameter values for which this phenomenon arises are biologically realistic remains a conjecture as field data will be needed to parameterize such occurrence. The Centre Manifold theory was used to determine the local asymptotic stability of the endemic equilibrium. Our results provide a perspective for understanding the complexity of MDR-TB, and the model can be applied in most settings where MDR-TB is present.

Numerical simulations suggest that, in the absence of any intervention, both TB strains cannot be eliminated from the population as 𝑅0>1, and the disease persists at an endemic equilibrium. A critical factor in addressing MDR-TB is primary prevention through DOTS and management of patients requiring second-line drug-regimen. Treatment of latent and infectious individuals with sensitive TB showed that ordinary TB can be completely eradicated. Thus, effective treatment for latent and infectious individuals with ordinary TB results in a reduction of MDR-TB, since the emergence of most MDR-TB cases is due to failure to provide TB drugs on time, as identifying latently infected individuals with sensitive and putting them on treatment is crucial in reducing new cases of resistant TB [29]. Also effective chemoprophylaxis and treatment of infectives result in a reduction of MDR-TB cases since most MDR-TB cases are a result of inappropriate treatment [5]. Treatment for infectious individuals with MDR-TB alters TB epidemics because it reduces the spread of MDR-TB strains and this supports the analytical results. Hence, a decrease in MDR-TB cases implies a decrease in MDR-TB-related deaths as MDR-TB kills more people than ordinary TB.

Diagnosis also plays an important role in MDR-TB reduction. As the proportion of TB patients being presented for diagnosis is increased, the rate of treatment should be correlated to the number of diagnosed infected individuals so as to reduce the burden of TB [6]. Significant increase in the detection rate of infectious individuals in Nigeria has been recommended because DOTS failed to reduce the incidence rate in the country due to failure to adequately detect a huge number of active TB cases which are primarily responsible for the spread of the infection [30]. As more people go for TB diagnosis, MDR-TB decreases due to the fact that those diagnosed with the disease are placed on DOTS. Drug-resistant TB will remain a serious threat to our communities as long as many members of our society do not have regular access to medical care [31]. Health education is another important aspect in the fight against MDR-TB. Results showed that, if more people receive health education, then the burden of MDR-TB can be reduced since MDR-TB cases also arise due to noncompliance with TB treatment. Information/awareness campaigns are viable in order to sensitize people on the importance of completing their TB dosages. Despite the role of the control strategies in reducing the burden of MDR-TB, numerical simulations also show that at the current level of TB treatment, diagnosis and health education, MDR-TB can only be reduced significanly.

Incomplete data from the National TB Control Program proves to be a major challenge in deriving estimates for the key biological parameters to calibrate the model to Malawi. Nevertheless, resorting to the literature, fundamental parameters values mimicking the epidemic in the region were used as a basis for illustration. Although several assumptions are made in the process, our results are driven by the model formulation and its structure; however, they are applicable to the Malawi context and other settings with similar epidemic trend. In summary, adequate treatment of sensitive TB will result in a reduction of MDR-TB in Malawi as most MDR-TB cases come from failure to properly administer TB drugs. Furthermore, diagnosis and health education of infectives with sensitive TB is important in the reduction of new MDR-TB cases due to adherence and compliance to treatment. Scaling up diagnosis, treating, and TB education will help in reducing the burden of the disease. Treatment rate of infected individuals should be correlated to the number of diagnosed individuals, and policies should be put in place to minimize loss to follow up. MDR-TB eradication remains a challenge to National Tuberculosis Control Programs in most developing countries, hence strengthening control strategies is paramount to curtailing TB spread, especially as the incidence rate of MDR-TB seems to be on the increase.

Finally, we identify some limitations of this study. A more realistic perspective could have been achieved by including, vaccination of susceptible population, immigrants, and new born; efficacy of MDR-TB drugs and information campaigns; controlling the disease with a possible minimal cost and side effects using control theory; estimating the cost-effectiveness of these control measures. Most parameter values are obtained from different sources giving rise to parameter uncertainty regarding their exact value. Our results which are driven by the model structure and its formulation are sensitive to the choice of parameter values. However, it is worth stressing that the main goal of this work was to provide a theoretical framework where the emergence of drug-resistant and MDR-TB can be addressed using a dynamical model. We focused on the population-level dynamics and potential benefits associated with implementation of various control strategies. It is our hope that the theoretical results obtained from this study will stimulate further interest in developing more complex models, be it agent based or network.

Appendices

A. Computation of the Effective Reproduction Number 𝑅𝑒

Following [19], the associated matrices 𝑖 for new infections terms and 𝒱𝑖 for the remaining transition terms are, respectively, given by 𝑖=1𝜌𝑠𝜆𝑆+𝑅𝑠𝜌𝑠𝜆𝑠𝑆1𝜌𝑟𝜆𝑆+𝑅𝑟𝜌𝑟𝜆𝑟𝑆00𝒱,(A.1)𝑖=𝛾𝑠𝜆𝑠+𝜆𝑟𝐸𝑠+𝜙𝑠+𝑘1𝐸+𝜇𝑠𝑞𝜙𝑟𝐼𝑠𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠𝐼+𝜎+𝜇𝑠+𝜆𝑟𝐼𝑠𝛾𝑠𝜆𝑠+𝑘1𝐸𝑠𝑘2𝐸+𝜇𝑟+𝛾𝑟𝜆𝑟+𝜆𝑠𝐸𝑟(1(𝑝+𝑞))𝜙𝑟𝐼𝑠𝑑𝑟+𝜙+𝜑𝑟𝐼+𝜇𝑟𝐸𝑠+𝐼𝑠𝜆𝑟𝜆𝑠+𝑘2+𝛾𝑟𝜆𝑟𝐸𝑟𝜎𝐼𝑠𝜆𝜇𝑆+𝑠+𝜆𝑟𝜆𝑆Λ𝜇𝑅+𝑠+𝜆𝑟𝑅𝑎𝑝𝜙𝑟𝐼𝑠𝜙𝑠𝐸𝑠𝜙𝐼𝑟.(A.2) Evaluating the partial derivatives of (A.1) at 𝐸0 and bearing in mind that system (1) has four infected classes, namely, 𝐸𝑠,𝐼𝑠,𝐸𝑟, and 𝐼𝑟, we obtain 0𝐹=1𝜌𝑠𝑐𝛽𝑠000𝜌𝑠𝑐𝛽𝑠000001𝜌𝑟𝑐𝛽𝑟000𝜌𝑟𝑐𝛽𝑟=𝐹100𝐹2,(A.3) where 𝐹1=01𝜌𝑠𝑐𝛽𝑠0𝜌𝑠𝑐𝛽𝑠𝐹2=01𝜌𝑟𝑐𝛽𝑟0𝜌𝑟𝑐𝛽𝑟.(A.4) Similarly, partial differentiation of (A.2) with respect to 𝐸𝑠,𝐼𝑠,𝐸𝑟, and 𝐼𝑟 at 𝐸0 gives 𝜙𝑉=𝑠+𝑘1+𝜇𝑞𝜙𝑟00𝑘1𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇000(1(𝑝+𝑞))𝜙𝑟𝑘2+𝜇00𝜎𝑘2𝑑𝑟+𝜙+𝜑𝑟=𝑉+𝜇10𝑉3𝑉2,(A.5)

where 𝑉1=𝜙𝑠+𝑘1+𝜇𝑞𝜙𝑟𝑘1𝑑𝑠+𝑎𝜙𝑟+𝜑𝑠,𝑉+𝜎+𝜇2=𝑘2+𝜇0𝑘2𝑑𝑟+𝜙𝑟+𝜑𝑟,𝑉+𝜇3=0(1(𝑝+𝑞))𝜙𝑟.0𝜎(A.6) The effective reproduction number, denoted by 𝑅𝑒, is given by 𝑅𝑒=𝜌(𝐹𝑉1), where 𝜌 denotes the spectral radius (or the dominant eigenvalue of matrix 𝐹𝑉1). The dominant eigenvalues of matrix 𝐹𝑉1 are given by 𝑅𝑠=𝐹1𝑉11=𝑐𝛽𝑠𝜇𝜌𝑠+𝜙𝑠𝜌𝑠+𝑘1𝜙𝑠+𝑘1𝑑+𝜇𝑠+𝑎𝜙𝑟+𝜑𝑠+𝜎+𝜇𝑞𝜙𝑟𝑘1,𝑅𝑟=𝐹2𝑉21=𝑐𝛽𝑟𝑘2+𝜇𝜌𝑟𝑘2𝑑+𝜇𝑟+𝜙+𝜑𝑟.+𝜇(A.7) Therefore, 𝑅𝑒=𝜌(𝐹𝑉1)=max{𝑅𝑠,𝑅𝑟}, where 𝑅𝑠 and 𝑅𝑟 are, respectively, the reproduction numbers for drug sensitive TB strain only and MDR-TB strain only. 𝑅𝑒 measures the average number of new infections generated by a typical infectious individual in a community where intervention strategies are in place. Thus, in the absence of diagnosis, treatment and health education (i.e., 𝜙𝑠=𝜙𝑟=𝜙=𝑎=𝜎=0), (A.7) reduces to 𝑅0𝑠=𝑐𝛽𝑠𝑘1+𝜇𝜌𝑠𝑘1𝑑+𝜇𝑠+𝜑𝑠,𝑅+𝜇0𝑟=𝑐𝛽𝑟𝑘2+𝜇𝜌𝑟𝑘2𝑑+𝜇𝑟+𝜑𝑟,+𝜇(A.8) and 𝑅0=max{𝑅0𝑠,𝑅0𝑟}.

B. Proof of the Stability of the EE

B.1. The Bifurcation Theorem

This Theorem is proven in Castillo-Chavez and Song [32].

Theorem B.1. Consider the following general system of ordinary differential equations with a parameter 𝜙: 𝑑𝑥𝑑𝑡=𝑓(𝑥,𝜙),𝑓𝑛×𝑛,𝑓2(𝑛×),(B.1) where 0 is an equilibrium point of the system, that is, 𝑓(0,𝜙)0forall𝜙 and(i)𝐴=𝐷𝑥𝑓(0,0)=[(𝜕𝑓𝑖/𝜕𝑥𝑗)(0,0)] is the linearization matrix of the system around the equilibrium 0 with 𝜙 evaluated at 0;(ii)zero is a simple eigenvalue of 𝐴 and all other eigenvalues of 𝐴 have negative real parts;(iii)matrix 𝐴 has a right eigenvector 𝑤 and a left eigenvector 𝑣 corresponding to the zero eigenvalue.

Let 𝑓𝑘 be the 𝑘th component of 𝑓 and𝑎=𝑛𝑘,𝑖,𝑗=1𝑣𝑘𝑤𝑖𝑤𝑗𝜕2𝑓𝑘𝜕𝑥𝑖𝜕𝑥𝑗(0,0),𝑏=𝑛𝑘,𝑖=1𝑣𝑘𝑤𝑖𝜕2𝑓𝑘𝜕𝑥𝑖𝜕𝜙(0,0).(B.2) The local dynamics of system (B.1) around 0 is governed by the signs of a and b.(i)𝑎>0, 𝑏>0. When 𝜙<0 with |𝜙|1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium; when 0<𝜙1, 0 is unstable and there exists a negative and locally asymptotically stable equilibrium.(ii)𝑎<0, 𝑏<0. When 𝜙<0 with |𝜙|1, 0 is unstable; when 0<𝜙1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium.(iii)𝑎>0, 𝑏<0. When 𝜙<0 with |𝜙|1, 0 is unstable, and there exists a locally asymptotically stable negative equilibrium; when 0<𝜙1, 0 is stable, and a positive unstable equilibrium appears.(iv)𝑎<0, 𝑏>0. When 𝜙 changes from negative to positive, 0 changes its stability from stable to unstable. Correspondingly a negative unstable equilibrium becomes positive and locally asymptotically stable.

Particularly if 𝑎>0 and 𝑏>0, then a backward bifurcation occurs at 𝜙=0.

B.2. Proof of Local Asymptotic Stability of 𝐸2

Again, using the Center Manifold theory [22], the local asymptotic stability of 𝐸2 is established. To this effect, the following change of variables is made; 𝑆=𝑥1,𝐸𝑟=𝑥2,𝐼𝑟=𝑥3,𝑅=𝑥4 (note that 𝐸𝑠=𝐼𝑠=0 at this point) so that 𝑁=𝑥1+𝑥2+𝑥3+𝑥4. Using vector notation 𝑋=[𝑥1,𝑥2,𝑥3,𝑥4]𝑇, the system (13) can be written in the form 𝑑𝑋𝑓𝑑𝑡=𝐹=1,𝑓2,𝑓3,𝑓4,(B.3) such that 𝑥1(𝑡)=𝑓1=Λ𝑐𝛽𝑟𝑥3𝑥1𝑁𝜇𝑥1,𝑥2(𝑡)=𝑓2=𝑐𝛽𝑟𝑥31𝜌𝑟𝑥1+𝑥4𝑁𝛾𝑟𝑐𝛽𝑟𝑥3𝑥2𝑁𝑘2𝑥+𝜇2,𝑥3(𝑡)=𝑓3=𝑐𝛽𝑟𝜌𝑟𝑥3𝑥1𝑁+𝛾𝑟𝑐𝛽𝑟𝑥3𝑥2𝑁+𝑘2𝑥2𝑑𝑟+𝜙+𝜑𝑟𝑥+𝜇3,𝑥4(𝑡)=𝑓4=𝜙+𝜑𝑟𝑥3𝑐𝛽𝑟𝑥3𝑥4𝑁𝜇𝑥4.(B.4) The Jacobian matrix of (B.4) at 𝐸0𝑟 is given by 𝐽𝐸0𝑟=𝜇0𝑐𝛽𝑟0𝑘02+𝜇1𝜌𝑟𝑐𝛽𝑟0𝑘2𝑐𝛽𝑟𝜌𝑟𝑑𝑟+𝜙+𝜑𝑟0+𝜇00𝜙+𝜑𝑟𝜇.(B.5) From (B.5), it can be shown that the the reproduction number is 𝑅𝑟=𝑐𝛽𝑟𝑘2+𝜇𝜌𝑟𝑘2𝑑+𝜇𝑟+𝜙+𝜑𝑟+𝜇.(B.6) If 𝛽𝑟 is the bifurcation point and if we consider the case when 𝑅𝑟=1 and then solve for 𝛽𝑟, we obtain 𝛽𝑟=𝛽=𝑘2𝑑+𝜇𝑟+𝜙+𝜑𝑟+𝜇𝑐𝑘2+𝜇𝜌𝑟.(B.7) System (B.5) with 𝛽𝑟=𝛽 has a simple zero eigenvalue, hence we can use the center manifold theory in the analysis of the dynamics of system (B.5) near 𝛽𝑟=𝛽. The Jacobian matrix (B.5) near 𝛽𝑟=𝛽 has a right eigenvector associated with the zero eigenvalue given by 𝑤=[𝑤1,𝑤2,𝑤3,𝑤4]𝑇, where 𝑤1=𝑐𝛽𝑤3𝜇,𝑤2=1𝜌𝑟𝑐𝛽𝑤3𝑘2=𝑑+𝜇𝑟+𝜙+𝜑𝑟+𝜇𝑐𝛽𝜌𝑟𝑤3𝑘2,𝑤3=𝑤3>0,𝑤4=𝜙+𝜑𝑟𝑤3𝑘2.(B.8)

The left eigenvector of (B.5) associated with the zero eigenvalue at 𝛽𝑟=𝛽 is given by 𝑣=[𝑣1,𝑣2,𝑣3,𝑣4]𝑇, where 𝑣1=𝑣4=0,𝑣3=𝑣3𝑣>0,2=𝑘2𝑣3𝑘2=𝑑+𝜇𝑟+𝜙+𝜑𝑟+𝜇𝑐𝛽𝜌𝑟𝑣31𝜌𝑟𝑐𝛽.(B.9)

We use Theorem  2.5 in [32] to find the conditions for the occurrence of backward bifurcation.

Computations of a and b
For system (B.4), the partial derivatives of 𝐹 associated with 𝑎 at 𝐸𝑟0 are given by 𝜕2𝑓2𝜕𝑥3𝜕𝑥2=𝜇𝑐𝛽1𝜌𝑟+𝛾𝑟Λ,𝜕2𝑓2𝜕𝑥3𝜕𝑥32=1𝜌𝑟𝜇𝑐𝛽Λ,𝜕2𝑓3𝜕𝑥3𝜕𝑥2=𝜇𝑐𝛽𝜌𝑟𝛾𝑟Λ,𝜕2𝑓3𝜕𝑥3𝜕𝑥3=2𝜌𝑟𝜇𝑐𝛽Λ.(B.10) Following (B.10), we have 𝑎=𝑣2𝑛𝑖,𝑗=1𝑤𝑖𝑤𝑗𝜕2𝑓2𝜕𝑥𝑖𝜕𝑥𝑗(0,0)+𝑣3𝑛𝑖,𝑗=1𝑤𝑖𝑤𝑗𝜕2𝑓3𝜕𝑥𝑖𝜕𝑥𝑗(0,0),=𝑣2𝑤3𝑤2𝜇𝑐𝛽1𝜌𝑟+𝛾𝑟Λ2𝑤231𝜌𝑟𝜇𝑐𝛽Λ+𝑣3𝑤3𝑤2𝜇𝑐𝛽𝜌𝑟𝛾𝑟Λ2𝑤23𝜌𝑟𝜇𝑐𝛽Λ,=2𝜇𝑐𝛽𝑣3𝑤23Λ𝑘2+𝜇2𝛾𝑟1𝜌𝑟𝜇𝑐𝛽1𝜌𝑟𝜇𝑐𝛽𝜌𝑟𝑘2+𝜇2𝜌𝑟𝑘2++𝜇1𝜌𝑟𝑘2.(B.11) We see from (B.11) that 𝑎<0 whenever 𝑚<𝑛 and 𝑎>0 whenever 𝑚>𝑛, where 𝑚=𝛾𝑟1𝜌𝑟𝜇𝑐𝛽,𝑛=1𝜌𝑟𝜇𝑐𝛽𝜌𝑟+𝑘2+𝜇2𝜌𝑟𝑘2++𝜇1𝜌𝑟𝑘2.(B.12) The nonzero partial derivatives of 𝐹 associated with 𝑏 at 𝐸𝑟0 are given by 𝜕2𝑓2𝜕𝑥3𝜕𝛽=1𝜌𝑟𝜕𝑐,2𝑓3𝜕𝑥3𝜕𝛽=𝑐𝜌𝑟.(B.13) It follows from (B.13) that, 𝑏=𝑣2𝑛𝑖=1𝑤𝑖𝜕2𝑓2𝜕𝑥𝑖𝜕𝛽(0,0)+𝑣3𝑛𝑖=1𝑤𝑖𝜕2𝑓3𝜕𝑥𝑖𝜕𝛽(0,0),=𝑣2𝑤31𝜌𝑟𝑐+𝑣3𝑤3𝜌𝑟=𝑐,𝑐𝑣3𝑤3𝑘21𝜌𝑟+𝜌𝑟𝑘2+𝜇𝑘2+𝜇>0.(B.14) Therefore, 𝑏>0 and 𝑎<0 or 𝑎>0 depending on whether 𝑚<𝑛 or 𝑚>𝑛. We have therefore, established the following result.

Theorem B.2. If 𝑚>𝑛, 𝑎>0, then model system (13) has a backward bifurcation at 𝑅𝑟=1, otherwise 𝑎<0 and a unique endemic equilibrium 𝐸2 is locally asymptotically stable for 𝑅𝑟>1 but close to 1.

B.3. Existence of Backward Bifurcation of the Full Model

From model system (1), we make the following change of variables, that is, 𝑆=𝑥1,𝐸𝑠=𝑥2,𝐼𝑠=𝑥3,𝐸𝑟=𝑥4,𝐼𝑟=𝑥5,𝑅=𝑥6, such that 𝑁=𝑥1+𝑥2+𝑥3+𝑥4+𝑥5+𝑥6. Further, by using vector notation 𝑋=[𝑥1,𝑥2,𝑥3,𝑥4,𝑥5,𝑥6]𝑇, system (1) can be written in the form 𝑑𝑋/𝑑𝑡=𝐹(𝑋), where 𝐹=(𝑓1,𝑓2,𝑓3,𝑓4,𝑓5,𝑓6) as follows: 𝑆(𝑡)=𝑓1𝜆=Λ𝑠+𝜆𝑟𝑥1𝜇𝑥1,𝐸𝑠(𝑡)=𝑓2=1𝜌𝑠𝑥1+𝑥6𝜆𝑠𝛾𝑠𝜆𝑠+𝜆𝑟𝑥2𝜙𝑠+𝑘1𝑥+𝜇2+𝑞𝜙𝑟𝑥3,𝐼𝑠(𝑡)=𝑓3=𝜌𝑠𝜆𝑠𝑥1+𝛾𝑠𝜆𝑠+𝑘1𝑥2𝜆𝑟𝑥3𝑑𝑠+𝑎𝜙𝑟+𝜎+𝜑𝑠𝑥+𝜇3,𝐸𝑟(𝑡)=𝑓4=1𝜌𝑟𝑥1+𝑥6𝜆𝑟+(1(𝑝+𝑞))𝜙𝑟𝑥3𝜆𝑠+𝛾𝑟𝜆𝑟+𝑘2𝑥+𝜇4,𝐼𝑟(𝑡)=𝑓5=𝜌𝑟𝜆𝑟𝑥1+𝜆𝑠+𝛾𝑟𝜆𝑟+𝑘2𝑥4+𝜆𝑟𝑥2+𝜆𝑟𝑥3+𝜎𝑥3𝑑𝑟+𝜙+𝜑𝑟𝑥+𝜇5,𝑅(𝑡)=𝑓6=𝑎𝑝𝜙𝑟+𝜑𝑠𝑥3+𝜙𝑠𝑥2+𝜙+𝜑𝑟𝑥5𝜆𝑠+𝜆𝑟𝑥6𝜇𝑥6,(B.15) where 𝜆𝑠=𝑐𝛽𝑠𝑥3/𝑁 and 𝜆𝑟=𝑐𝛽𝑟𝑥5/𝑁. The Jacobian matrix of system (B.15) at 𝐸0 is given by𝐽𝐸0=𝜙𝜇000000𝑠+𝑘1+𝜇1𝜌𝑠𝑐𝛽𝑠+𝑞𝜙𝑟0000𝑘1𝑔100000𝑔2𝑘2𝑔+𝜇3000𝜎𝑘2𝑔400𝜙𝑠𝑎𝑝𝜙𝑟+𝜑𝑠0𝜙+𝜑𝑟𝜇(B.16)

where 𝑔1𝑑=𝑠+𝑎𝜙𝑟+𝜎+𝜑𝑠+𝜇+𝑐𝛽𝑠𝜌𝑠,𝑔2=(1(𝑝+𝑞))𝜙𝑟,𝑔3=1𝜌𝑟𝑐𝛽𝑟,𝑔4𝑑=𝑟+𝜙+𝜑𝑟+𝜇+𝜌𝑟𝑐𝛽𝑟.(B.17) It can be shown that the eigenvalues of (B.16) are expressed in terms of 𝑅𝑒=max{𝑅𝑠,𝑅𝑟}, where 𝑅𝑠 and 𝑅𝑟 are the reproduction numbers of drug sensitive and drug-resistant TB only sub-models respectively as seen earlier. 𝑅𝑒=max{𝑅𝑠,𝑅𝑟} implies that the two TB strains (drug sensitive and drug-resistant) escalate each other. Thus, when the two reproduction numbers exceed unity, that is, 𝑅𝑠>1 and 𝑅𝑟>1, there is always coexistence (endemic case) of these two strains regardless of which reproduction number is greater as shown in Theorem 4. If 𝑅𝑒=max{𝑅𝑠,𝑅𝑟}, then from Theorem 2, the drug sensitive TB only sub-model has a backward bifurcation for values of 𝑅𝑠 such that 𝑅𝑐𝑠<𝑅𝑠<1 and Theorem 3 showed that the drug-resistant TB only sub-model exhibits backward bifurcation for values of 𝑅𝑟 such that 𝑅𝑐𝑟<𝑅𝑟<1. Thus, the coexistence model of TB will also exhibit the phenomenon of backward bifurcation whenever 𝑅𝑒=1

Acknowledgments

M. Maliyoni and P. M. M. Mwamtobe acknowledge with thanks the financial support by the Norad’s Masters programme in Mathematical Modelling (NOMA) at the University of Dar es Salaam, Tanzania, and the University of Malawi, The Polytechnic, for study leave. The authors thank the reviewers for insightful comments.