Mathematical analysis of a two-strain tuberculosis model in Bangladesh

Tuberculosis (TB) is an airborne infectious disease that causes millions of deaths worldwide each year (1.2 million people died in 2019). Alarmingly, several strains of the causative agent, Mycobacterium tuberculosis (MTB)—including drug-susceptible (DS) and drug-resistant (DR) variants—already circulate throughout most developing and developed countries, particularly in Bangladesh, with totally drug-resistant strains starting to emerge. In this study we develop a two-strain DS and DR TB transmission model and perform an analysis of the system properties and solutions. Both analytical and numerical results show that the prevalence of drug-resistant infection increases with an increasing drug use through amplification. Both analytic results and numerical simulations suggest that if the basic reproduction numbers of both DS (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{R}}_{{0{\text{s}}}}$$\end{document}R0s) and DR (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{R}}_{{0{\text{r}}}}$$\end{document}R0r) TB are less than one, i.e. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\max \left[ {{\text{R}}_{{0{\text{s}}}} ,{\text{ R}}_{{0{\text{r}}}} } \right] < 1,$$\end{document}maxR0s,R0r<1, the disease-free equilibrium is asymptotically stable, meaning that the disease naturally dies out. Furthermore, if \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{R}}_{{0{\text{r}}}} > {\text{max}}\left[ {{\text{R}}_{{0{\text{s}}}} ,1} \right]$$\end{document}R0r>maxR0s,1, then DS TB dies out but DR TB persists in the population, and if \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{R}}_{{0{\text{s}}}} > {\text{max}}\left[ {{\text{R}}_{{0{\text{r}}}} ,1} \right]$$\end{document}R0s>maxR0r,1 both DS TB and DR TB persist in the population. Further, sensitivity analysis of the model parameters found that the transmission rate of both strains had the greatest influence on DS and DR TB prevalence. We also investigated the effect of treatment rates and amplification on both DS and DR TB prevalence; results indicate that inadequate or inappropriate treatment makes co-existence more likely and increases the relative abundance of DR TB infections.

www.nature.com/scientificreports/ prevalence and incidence in India and the United States of America (USA) 13 . Kim et al. (2014) developed a mathematical model for TB with exogenous reinfection, examining the present situation of active TB incidence in Korea 14 . Liu et al. (2010) developed a TB model with seasonality to describe TB incidence rates with periodic properties in a mainland city of China 15 . A 10-compartment TB model constructed by Trauer et al. (2014) modelled limited vaccine effectiveness, reinfection, DR TB, and de novo resistance through treatment 16 . This study showed that the model could not be calibrated to the projected incidence rate without allowing for reinfection, which was modelled as a reversion to early latency, which has a higher rate of progression to disease compared with late latency. A theoretical framework for TB control strategies is developed by Blower et al. (1996) to assess the effect of chemoprophylaxis and treatment for eradication of TB epidemics. Results show that increasing chemoprophylaxis and treatment rates will decrease the severity of TB epidemics and reduce the probability that a latently infected individual will progress to active TB. This study recommended that in developing countries treatment failure rates must be lower than in developed countries in order to control TB outbreaks 17 . Kuddus et al. (2020) developed a two-strain TB model to determine the optimal time-varying combination of distancing, latent case finding, case holding and active case finding control strategies that minimize tuberculosis incidence over a fixed timespan. The major finding of this study was that combining two or more intervention strategies is more cost-effective compared to single-intervention strategies 12 . In 2010, Huo et al. developed a two-strain transmission dynamic model with susceptible, exposed and infected compartments that included drug-sensitive and drug-resistant TB 18 . Dynamical system analysis was used to ascertain the global stability of the disease-free equilibrium and a mono-existent equilibrium. Results show that the drug-sensitive and drug-resistant TB will die out if the both basic reproduction numbers are less than one, otherwise an epidemic occurs 18 .
TB is an ancient disease, which is still a major public health problem in Bangladesh. The problem is aggravated by the increasing drug-resistant mutations that arise through inadequate or poor quality treatment. To evaluate the threat of genetic variations of DS and DR TB strains, we present a two-strain (DS and DR TB) Susceptible-Latent-Infectious-Removed-Susceptible (SLIRS) TB model and examine the development and transmission of DR TB. We acknowledge the possibility that an individual's position swaps from DS TB at primary appearance to DR TB at follow-up. This is the manner by which DR TB typically emerges in a population and is intended to replicate the phenotypic phenomenon of acquired drug resistance, known as amplification. The model can be used to examine the co-existent or competitive phenomena between DS and DR TB strains.
In this study, we perform a rigorous analytical and numerical analysis of our novel two-strain TB model properties and solutions from both the mathematical and biological perspectives. For an individual, we apply the next generation matrix technique (which is effectively a table describing the number of new infections generated by each individual infected with a particular strain) to determine analytic expressions for the basic reproduction numbers of the DS TB and DR TB strains. We find that these are key determinates for regulating model dynamics. We also present the mandatory conditions for the stability of the infection-free state, as well as infection mono-existence and co-existence.
To enhance and verify the analytical analysis, we apply numerical procedures to solve the system equations and investigate the dynamic epidemic trajectory for a variety of initial conditions and plausible parameter values. From the analytical and numerical viewpoints, the local and global stability of the infection-free equilibrium and mono-existent equilibrium are examined through standard dynamical systems analysis techniques 19 . The coexistent endemic equilibrium is also examined numerically. Following this, we implement a sensitivity analysis to examine the model parameters that have the greatest effect on DS, DR and total TB prevalence.

Methods and materials
Model formulation. We developed a deterministic mathematical model of the transmission of DS and DR TB strains between the following mutually exclusive compartments: susceptible S(t) , uninfected individuals who are susceptible to TB infection; those exposed to TB that become latently infected L i (t) (where the subscript i = s, r refers to quantities associated with drug-susceptible, s , and drug-resistant, r , TB infection), representing those who are infected and have not yet developed active TB; infectives I i (t) , comprising individuals with active TB that are capable of transmitting the infection; the recovered R(t) , who were previously infected and successfully treated and assumed to be temporarily immune to reinfection. We assume that DR TB is initially generated through inadequate and poor treatment of DS TB and could subsequently be transmitted to other individuals. Individuals may also return to the susceptible compartment following recovery at the constant per-capita rate γ due to loss of immunity. The total population size N(t) , is given by Individuals in the different compartments suffer from natural death at the same constant rate µ and active TB cases in I i (i = s, r) experience disease-related death at a rate φ i (i = s, r) . To ensure the population size remains constant, we replace all deaths as newborns in the susceptible compartment. Individuals in the S compartment may be infected with a circulating MTB strain i (i = s, r) at a time dependent rate i (t) = β i I i (t) where β i is the transmission rate between infected and susceptible individuals. Once infected with strain I i (i = s, r) , individuals move to the latently infected compartment L i . A proportion of those with latent infections progress to active TB as a result of endogenous reactivation of the latent bacilli at rate α i . Individuals with drug sensitive and DR active TB, I i , may eventually be detected and treated at rates τ s and τ r for DS and DR TB, respectively. A proportion (1 − ρ) of the treated DS active TB recover to move into the recovered compartment R , and the complementary proportion ρ develop drug resistance due to incomplete treatment or lack of strict compliance in the use of firstline drugs (drugs used to treat the DS forms of TB) to move into compartment I r . Furthermore, individuals  www.nature.com/scientificreports/ recover naturally at a rate ω i , moving from I i to R . The model flow diagram is presented in Fig. 1. The parameters and their values are given in Table 1. The transmission of DS and DR TB is given by the following deterministic system of nonlinear ordinary differential equations that describe the model: with non-negative initial conditions for the system of differential equations above, it is easy to express that each of the state variables remain non-negative for all t > 0 . Moreover, from Eqs. (2)-(7), we find that the total population size, N(t) satisfies Given the constant population size and positivity of solutions it generally follows that each of the compartment states S, L, I, etc. are bounded. Therefore, based on the considerations, the feasible region for the system Eqs.
(2)-(7) is given by Estimation of model parameters. Most parameters were estimated from the global literature or taken from National TB Control Program (NTP) evaluation (Table 1). We determined the total population of Bangladesh in 2015 is about 159,000,000 based on the National TB Control Program (NTP) reports. The natural death rate (µ) was considered as the inverse of life expectancy (70 years) of Bangladesh. According to the data from the WHO report 20 the treatment success rates for DS (τ s ) and DR (τ m ) TB in Bangladesh were about 94% and 78% respectively. We estimated two critical residual model parameters, β s and β m by fitting our model to the DS and DR TB incidence data in Bangladesh taken from the WHO and NTP report from 2000 to 2018 1,5,21 . The best fitted (see Fig. 2) parameter values β s and β m were obtained, using the least-square optimization method to dI r dt = α r L r − ω r I r − µI r + ρτ s I s − φ r I r − τ r I r , (c) D = (S, L s , I s , L r , I r , R) ∈ R 6 + : S + L s + I s + L r + I r + R = N . Sensitivity analysis. We calculated the Partial Rank Correlation Coefficients (PRCCs) between each of the model parameters and several outcome measures of significance, a global sensitivity analysis method using Latin Hypercube Sampling (LHS). Explicitly, a uniform distribution is allocated from half to fourfold the baseline value (see Table 1) for each model parameter (i.e.β s , α s , ω s , φ s , τ s , β r , α r , ω r , φ r , τ r , ρ and γ ) and a total of 10,000,000 random draws are taken for each. The model is then simulated for each of the 10,000,000 parameter sets and relevant outputs such as disease prevalence and incidence are recorded. Here the model outputs we considered were DS TB (I s ) , DR TB (I r ) and total TB (I s + I r ) prevalence at equilibrium. From the analytical formula for R 0s and R 0r , the sensitivity indices (measuring the relative change in a variable when a parameter changes) ϒ i j can be derived applying the technique in 22 to each of the model parameters. For example, for β s we have: Here ϒ R 0s β s is the sensitivity index for the basic reproduction number R 0s , as we vary β s .
Ethical approval. The study was based on aggregated DS and DR TB surveillance data taken from the World Health Organization and the National TB Control Program in Bangladesh 1,5,21 . No confidential information was involved because analyses were performed at the aggregate level. All of the methods were conducted in accordance with the approved research protocol. The research protocol was approved by the James Cook University human ethics approval board, H7300.

Results
Basic reproduction number. The model has four infected states: L s , I s , L r , I r , and two uninfected states: S and R. At the infection-free steady state, To calculate the basic reproduction numbers of the DS and DR TB strains we follow 23 and focus on the linearized infection subsystem derived from Eqs. (2)- (7): where χ s = ω s + φ s + τ s + µ and χ m = ω r + φ r + τ r + µ are the total removal rates from the DS and DR active TB infection states respectively.
Here, the ODEs (8)- (11) represent the changes in the infected states about the infection-free equilibrium (that is, when assuming the reduction in the susceptible population as a result of infection is negligible).
By setting X T = (L s , I s , L r , I r ) T , where T denotes transpose, now we can write the infection subsystem in the following form The matrix T corresponds to transmissions and the matrix to transitions. They are acquired from system (8) − (11) by splitting the transmission events from other events. If we indicate to the infected states with indices i and j, with i, j ∈ 1, 2, 3, 4, then entry T ij is the rate per unit time at which persons in infected state j give rise to persons in infected state i in the system. Hence, for the subsystem (8) − (11) we obtain The next generation matrix, K , is given by 24 (note the essential minus sign) dI r dt = α r L r − χ r I r + ρτ s I s , www.nature.com/scientificreports/ The leading eigenvalues of the next generation matrix (K) are the basic reproduction numbers for DS and DR TB; they describe the average number of secondary infections generated by a single infected individual. Therefore the basic reproduction number for DS and DR TB are: and Here α s (αs+µ) and α r (αr+µ) are the probability of transitioning from the latent compartment to the infectious compartment of the DS and DR strains respectively. Further, 1 χ s and 1 χ r represent the time spent by infectious individuals in states I s and I r respectively.
The quantities R 0s and R 0m represents the expected number of secondary cases for DS and DR TB produced by the single infectious introduced into a completely susceptible population. Remarkably we identify that the basic reproduction numbers of DS TB (R 0s ) and DR TB (R 0r ) are completely independent of the proportion of amplification ρ 24 .
We provide detailed analysis of the proposed TB model (2)- (7) for the existence and stability of disease-free (E * ) , mono-existent (E ∧ ) and co-existent E † endemic equilibrium points in the supplementary materials: existence of equilibria section and stability analysis section.
In summary the disease-free equilibrium is globally asymptotically stable when max[R 0s , R 0r ] < 1 . We demonstrate that the mono-existent equilibrium is locally stable provided R 0r > max[1, R 0s ] (see supplementary materials: stability analysis section).
To establish the nature of the co-existence equilibrium, we implemented numerical analysis using the Monte Carlo method in 25 to verify the conditions R 0s > max[R 0r , 1] by calculating the real part of the eigenvalues of the Jacobian matrix evaluated at the co-existent endemic equilibrium E † . According to E † the coordinates of the steady states can be expressed in terms of twelve of the model's parameters (α s , α r , β s , β r , ω s , ω r , φ s , φ r , τ s , τ r , ρ, γ) whose baseline values are given in Table 1. The sampling pool Q =  Table 1) and θ > 0 defines the variation range. Our sampling comprised 10,000 confounding scenarios Q = q 1 , . . . , q 12 ∈ Q where each q i ∈ M i , i = 1, . . . , 12 was randomly chosen for θ = 0.20 (i.e. 20% deviation from the baseline values) under a uniform distribution with no correlation between parameters. The simulation results are presented in Fig. 3. Figure 3 showed that the co-existent endemic equilibrium is stable (yellow-coloured points) when R 0s > max[R 0r , 1] , otherwise unstable (blue-coloured points). Figure 4 represents the existence and local stability regions for the disease-free (E * ) , mono-existent (E^) and co-existent E † endemic equilibrium points corresponding to the drug-susceptible and drug-resistant TB basic reproduction numbers R 0s and R 0r . The magenta shaded region indicates the disease-free equilibrium which is bounded by max[R 0s , R 0r ] < 1 . The green shaded area illustrates the mono-existent equilibrium where R 0r > max[1, R 0s ] . Finally, the co-existent equilibrium in the yellow shaded region where R 0s > max[R 0r , 1].
Various numerical simulations were carried out using the Matlab programming language to support the analytic results and to observe the impact of amplification and the DS TB treatment rate on total and DR TB prevalence. We applied different initial conditions for both DS and DR TB of all populations and obtained the stability results for the model equilibria, finding TB disease will ultimately die out from the population when the condition max[R 0s , R 0r ] < 1, holds. The condition R 0r > max[R 0s , 1] implies that DS TB dies out but DR TB persists in the population. Furthermore, the condition R 0s > max[R 0r , 1] suggests that both DS TB and DR TB persist in the population (see supplementary materials: stability analysis section, Figures S1, S2 and S3).
Simulating the impact of public health effects through changes in model parameters. Public health interventions can be modelled by changing key parameters. For example, reducing contact between people can be modelled by changing β s or β r , detecting and treating infectious people can be modelled by increasing τ s or τ r . The resulting reproduction number after a change in a public health intervention or after a substantial proportion of the population is removed from the susceptible pool is often called the effective reproduction number. Figure 5 shows the impact of amplification (ρ) on DS and DR TB prevalence in the first region ( ρ < ∼ 0.6 ), DS TB initially dominates, but for ρ > ∼ 0.6 the DR TB has higher prevalence due to the amplification pathway when R 0s > max[R 0r , 1] . These outcomes implicate the use of drug treatment in the context of DR TB and note that treatment of DS TB can provoke the emergence of DR TB, even if it has lower reproductive fitness than the susceptible strain, including when its reproduction number is less than one. Further, under these situations, our investigation reveals that this emergence of DR TB will be overcome if the treatment rate is adequate to eliminate DS TB from the population (i.e. is able to achieve an effective reproduction number of DS TB less than one). Figures 6 and 7 display the impact of the DS TB treatment rate (τ s ) and amplification ( ρ ) on the total TB and DR TB prevalence, when both infectious rates (β s , β r ) are fixed. When we increase the amplification, the total TB and DR TB prevalence also increases. Conversely, Fig. 7 shows that for high amplification, DR TB prevalence increased when the DS TB treatment rate increased from zero to around 0.8 to 0.9 and then dropped to a common point. For small amplification values, the DR TB proportion only raised up to the common point. This point is the DR TB-only equilibrium which becomes stable when the basic reproduction ratio of DR TB becomes higher than the effective reproduction number of DS TB and higher than one. Numerical analysis shows that (a) R 0s = Nα s β s (α s + µ)χ s , www.nature.com/scientificreports/ for adequately high amplification, the DR TB prevalence will exceed that of its inherent equilibrium value when the DS TB is in existence and is being treated. From the above numerical analysis, it is clear that proper treatment is very important for DS TB patients otherwise it will lead to the creation of new cases of DR TB which can become the dominant strain. We performed sensitivity analysis to explore the quantitative relationship between parameters and model outcomes including DS, DR and total TB prevalence. Figure 8 displays the association between the co-existent equilibrium value of DS TB prevalence ( I s ) and model parameters β s , α s , ω s , φ s , τ s , β r , α r , ω r , φ r , τ r , ρ and γ , when R 0s > max[R 0r , 1] . From Fig. 8 it is observed that DS TB prevalence ( I s ) has a positive association with β s , α s , β r , α r , and γ , indicating that a positive change in any of these parameters will increase the DS TB prevalence ( I s ). We found that the association between the transmission rate β s and DS TB prevalence ( I s ) has the highest positive impact of all parameters, which is expected given the reproduction number is proportional to β s . The public health lesson is that reducing exposure between infected and uninfected individuals is the most important public health intervention, dominating the importance of good treatment, reduced relapse, and amplification. Public health services should concentrate on those who are in close and prolonged contact with an infectious individual (i.e. individuals who are living in the same household with infectious TB cases) who currently have extensive exposure.
In contrast, parameters ω s , φ s , τ s , ω r , φ r , τ r and ρ have a negative correlation with DS TB prevalence ( I s ), which means that increasing these parameters values will consequently decrease the I s prevalence. From the analysis, we also found that the treatment rate (τ s ) is the second most important parameter and its association with DS TB prevalence. Our finding is consistent with previous studies 12,14,26 , and the WHO which recommended that early treatment and cure of infectious cases of TB to cut the chain of transmission of TB infection in the community. Therefore, quick identification of presumptive TB cases, rapid diagnosis, and early initiation of treatment and successful completion of treatment are the most effective ways of preventing TB. Further, amplification has a negative impact on DS TB prevalence because some of the DS TB-infected individuals move to the DR TBinfected state due to incorrect treatment. Figure 9 represents the correlation between the DR TB prevalence and corresponding model parameters β s , α s , ω s , φ s , τ s , β r , α r , ω r , φ r , τ r , ρ and γ when R 0s > max[R 0r , 1] . Parameters β s , α s , β r , α r , ρ and γ have positive PRCC values and parameters ω s , φ s , τ s , ω r , φ r and τ r have negative PRCC values. Our finding is in line with observations and shows that amplification (ρ) has a positive correlation with DR TB prevalence at the co-existent equilibrium because it is the pathway from DS TB to DR TB. Figure 3. Stability of the co-existent endemic equilibrium as a function of R 0s /R 0m (as determined by the eigenvalues of the system Jacobian matrix). Random parameter draws leading to eigenvalues with exclusively negative real parts (i.e. stable) are colored yellow, whilst those leading to eigenvalues with at least one positive real part (i.e. unstable) are colored blue. All stable points (yellow) lie above the line R 0s = R 0r , and all unstable points (blue) lie below.    www.nature.com/scientificreports/ We further performed sensitivity analysis for the total TB prevalence (I s + I r ) and DR TB prevalence according to their corresponding model parameters when R 0s > max[R 0r , 1] and R 0r > max[R 0s , 1] respectively (see supplementary materials: sensitivity analysis section, Figure S4 and S5).
The sensitivity indices of R 0s and R 0r with respect to each parameter are given in Table 2. In the sensitivity indices of R 0s and R 0r , the most significant parameters are the transmission rates of DS TB,β s and DR TB,β r . Since ϒ R 0s β s = 1 , and ϒ R 0r β r = 1 , increasing (or decreasing) the transmission rates, β s and β r of DS TB and DR TB by a particular percentage, increases (or decreases) the reproduction numbers R 0s and R 0r by the same percentage.

Discussion and conclusion
In this paper, we formulated and analysed a novel two-strain TB model with amplification: one strain for DS TB; and another for DR TB. Here, we considered amplification as the process by which an individual infected with DS TB develops infection with a resistant strain of TB, reflecting treatment failure for individuals on first line drug therapy.
We found three equilibrium points of our proposed model; the disease-free equilibrium; the mono-existent equilibrium, when DR TB dominated in this system; and the co-existent equilibrium, when DS TB dominated in this system. The next generation matrix method was used to calculate the basic reproduction number of the different TB strains, denoted by R 0s for DS TB and R 0r for DR TB. The value of the basic reproduction numbers, namely R 0s and R 0r , and biological parameters of the model, were estimated on the basis of available data and are tabulated in Table 1. Furthermore, dynamical system analyses were also used to investigate the local stability of the disease-free equilibrium and mono-existent equilibrium. This analysis showed that stability depends on the threshold quantities, i.e. the basic reproduction numbers R 0s and R 0r . If max [R 0s , R 0r ] < 1 , the disease-free equilibrium is globally asymptotically stable, which means that the disease naturally dies out. If R 0r > max[R 0s , 1] , DS TB dies out but DR TB persists in the population. If R 0s > max[R 0r , 1] , then DS TB and DR TB both persist in the population. This analysis can help us to identify regions in the parameter space where the various asymptotic states are stable or unstable, thus allowing us to predict the long-term behaviour of the DS and DR TB dynamics. This information can advise the Bangladesh National TB Control Program to reduce the period of infectiousness down until max [R 0s , R 0r ] < 1.
The existence and stability of the transmission dynamics of sensitive and resistant TB bacteria within and between-hosts models were considered in previous modelling studies 27,28 . Previous studies showed that if the basic reproduction numbers of sensitive bacteria and resistant bacteria is less than one then the bacteria are cleared. If the basic reproduction number of resistant bacteria is greater than one but that of the sensitive bacteria is less than one then only resistant bacteria persist. Finally, the previous studies also showed that both sensitive Figure 6. Effect of DS TB treatment rate ( τ s ) on total TB prevalence when β s and β r are fixed. All remaining parameters consider their baseline values as reported in Table 1 www.nature.com/scientificreports/ and resistant bacteria persist if the sensitive bacteria's basic reproduction number is greater than one and the resistant bacteria's basic reproduction number bigger than certain threshold but less than one 27,28 , which is not the same as our results, which do not require the basic reproduction number of the resistant bacteria to be above a threshold in order to achieve co-existance, owing to the amplification pathway. Several TB modelling studies 14,17 explored the impact of single and combination intervention strategies. Results show that combination intervention strategies are the most effective for reducing the burden of TB. In this study, we considered a Bangladesh-specific six-compartmental two-strain SLIRS model with amplification. An important feature of this model is the coupling between the two strains representing the flow of infected individuals who acquire resistance during treatment. We assumed amplification develops mainly through the choice of naturally occurring mutations in the presence of inappropriate treatment. We also considered a natural recovery parameter and treatment parameter which are not considered in the earlier analyses 14,17 . Here, we additionally allow recovered people to move to the susceptible class due to loss of immunity, a feature that is also not considered in previous studies 14,17,18 .
Our model determined that from the explicit formulae for R 0s and R 0r , it is clear that these basic reproduction numbers depend on transmission rates β s (β r ) , progression rates α s (α r ) , recovery rates ω s (ω r ) , disease related death rates φ s (φ r ) , and treatment rates τ s (τ r ) . From the sensitivity analysis it is also clear that the most important parameters are the transmission rates β s (β r ) followed by the treatment rates τ s (τ r ) . Therefore, to control and eradicate DS TB and DR TB infection, it is important to consider the following strategies: the first and most important strategy is to minimize the contact rates β s (β r ) with infected individuals. There are many ways we can minimize the contact rates, such as personal respiratory protection; this includes wearing of masks by patients to minimise dispersal of TB bacilli when they talk, cough, yawn or sneeze. In addition, basic infection control measures have to be taught to patients such as covering the nose and mouth during coughing and sneezing and Figure 9. PRCC values showing the association between model output I r and the model parameters β s , α s ω s , φ s , τ s , β r , α r , ω r , φ r , τ r andρ , when R 0s > R 0r and R 0s > 1. www.nature.com/scientificreports/ to discard used tissues into covered bins. Environmental control measures maximise dilution and air exchange and decontaminate air when adequate ventilation cannot be reached in high-risk areas. Any ventilation system must be monitored and maintained on a regular schedule. Diagnosis campaigns are also needed to control transmission. The second-most important strategy is to increase the treatment rates τ s (τ r ) of infective individuals. Moreover, higher values of amplification (ρ) when combined with high τ s , lead to higher prevalence of DR TB-simulating inadequate or poorly administered treatment resulting in the creation of new cases of DR TB. Therefore, we suggest that for problematic DR TB, measures of the reproduction numbers of the DS and DR TB be performed together with the risk of amplification, to assure optimal levels of treatment be utilized to reduce the risk of DR TB. Well-administered first-line treatment for DS TB is the best way to prevent acquisition of resistance. Timely identification of DR TB and adequate treatment regimens with second-line drugs administered are essential to prevent rises in DR TB prevalence. However, in developing countries (e.g. Bangladesh) it is very difficult to isolate infectious individuals due to the high cost of long-term treatment. Therefore, we propose the most feasible and optimal strategy to eliminate DS TB and DR TB in Bangladesh is to increase the treatment rates by decreasing the treatment cost allowing universal health care.