Modeling Transmission of Tuberculosis with MDR and Undetected Cases

Abstract

This paper presents a novel mathematical model with multidrug-resistant (MDR) and undetected TB cases. The theoretical analysis indicates that the disease-free equilibrium is globally asymptotically stable if 𝑅₀<1; otherwise, the system may exist a locally asymptotically stable endemic equilibrium. The model is also used to simulate and predict TB epidemic in Guangdong. The results imply that our model is in agreement with actual data and the undetected rate plays vital role in the TB trend. Our model also implies that TB cannot be eradicated from population if it continues to implement current TB control strategies.

1. Introduction

China is among the highest TB burden countries in the world and is second only to India. The number of MDR-TB cases is about one third of that of the world [1, 2]. For an MDR-TB patient, the cost for the treatment is usually 10–100 times higher than that of a common TB patient [3]. Moreover, undetected TB patients pose threat to others, which plays a major role in mycobacterium transmission among the general population [4]. It is no doubt that the earlier the diagnosis and treatment, the better in attenuating the transmissions. Clearly, the understanding of the mechanism and contribution of undetected patients in transmission is important to conduct targeted interventions and control TB decease.

Blower et al. proposed a transmission dynamics among the cases infected by drug-sensitive and the drug-resistant strains, respectively [5]. Many subsequent investigations focused on similar topics with various perspectives about the role of drug-resistant strain [3, 6–8]. However, these works failed to independently identify the role of the undetected patients in TB transmission. motivated by the previous work, we carefully divide the infectious cases into undetected cases and timely detected cases in this paper. This improved classification enables us to interpret the real situation in a better manner. To analyze the stability of the disease-free equilibrium in our transmission model, we introduce a novel Lyapunov function. We prove that the disease-free equilibrium is globally asymptotically stable when the basic reproduction number is less than one, and a forward (transcritical) bifurcation exists otherwise.

While Guangdong is a developed province of China, the number of TB cases in Guangdong is about 10% of that of national TB cases. Mass migrants, MDR, and HIV/AIDS coinfections make the situation of TB transmission very complex in this region. We apply our model to describe the TB situation in Guangdong and simulate the distribution and predication. The simulations indicate that undetected cases do play an important role in TB transmission.

The remainder of the paper is organized as follows. The proposed TB model is formulated in Section 2. Stability criteria of the equilibria are presented in Section 3. Data analysis and simulation studies are provided in Section 4. Finally, Section 5 concludes with discussions and conclusions.

2. Model Description

In this paper, we consider the transmission diagram shown in Figure 1 and the mathematical model is given by 𝑑𝑆𝑑𝑡=𝐴−𝛽₁𝐼_𝑈𝑆−𝛽₂𝐼_𝐺𝑆−𝛽₃𝐼_𝑅𝑆−𝜇𝑆,𝑑𝐸𝑑𝑡=𝛽₁𝐼_𝑈𝑆+𝛽₂𝐼_𝐺𝑆+𝛽₃𝐼_𝑅𝑆−𝜇+𝑘+𝑘₁𝐸+𝜔𝑅,𝑑𝐼_𝑈𝑑𝑡=𝑘𝐸−𝜇+𝜇_𝐷+𝜐𝐼_𝑈,𝑑𝐼_𝐺𝑑𝑡=𝜐𝑚𝐼_𝑈+𝑘₁𝑝𝐸−𝜇+𝜇^′_𝐷+𝑟₁𝐼_𝐺−𝜎𝐼_𝐺,𝑑𝐼_𝑅()𝐼𝑑𝑡=𝜐1−𝑚_𝑈+𝑘₁1−𝑝𝐸−𝜇+𝜇_𝐷^′′+𝑟₂𝐼_𝑅+𝜎𝐼_𝐺,𝑑𝑅𝑑𝑡=𝑟₁𝐼_𝐺+𝑟₂𝐼_𝑅−𝜇+𝜔𝑅, (2.1) where 𝑆, 𝐸, 𝐼_𝑈, 𝐼_𝐺, 𝐼_𝑅, and 𝑅 represent the numbers of cases of susceptible, latent, undetected infectious TB, detected drug-sensitive infectious TB, detected infectious MDR-TB, and the treated, respectively. 𝐴 is the recruitment rate. 𝛽₁, 𝛽₂, and 𝛽₃ are the transmission rates of the undetected infectious, detected infectious drug-sensitive, and detected infectious MDR-TB cases for susceptible individuals, respectively. 𝜇 is the per capita natural mortality rate. 𝜇_𝐷, 𝜇^′_𝐷, and 𝜇_𝐷^′′ are the death rates (associated to TB) of undetected infectious, detected infectious drug-sensitive, and detected infectious MDR-TB cases, respectively. (𝑘+𝑘₁) is progressive rate from latent to infectious, with 𝑘₁ being the detected rate from latent cases to undetected cases in two weeks (from onset of symptoms to treatment, see [2]). The infectivity of the infected patient is strongest in the first two weeks after being infected. 𝜐 is the per capita detected rate from undetected to detected, called the delayed detected rate. 𝑟₁, 𝑟₂ are the treatment rates of detected infectious drug-sensitive and detected infectious MDR-TB cases, respectively. 𝜔 is the reactivation rate from treated to latent. 𝜎 is the treatment failure rate of detected drug-sensitive infectious patients who have been incompletely treated. We assume that patients without incomplete treatment will develop into MDR-TB as in [5] and the infectious transmission rate is bilinear.

Figure 1 

Transfer diagram of the model.

Let 𝑁(𝑡)=𝑆+𝐸+𝐼_𝑈+𝐼_𝐺+𝐼_𝑅+𝑅 be the total population. It follows from the system equations that 𝑑𝑁(𝑡)𝑑𝑡≤𝐴−𝜇𝑁, (2.2) which implies that lim_𝑡→∞𝑁(𝑡)≤𝐴𝜇. (2.3) Therefore, we have an attractor of system (2.1); that is, 𝑃=𝑆,𝐸,𝐼_𝑈,𝐼_𝐺,𝐼_𝑅,𝑅∈𝑅⁶₊∶𝑆+𝐸+𝐼_𝑈+𝐼_𝐺+𝐼_𝑅𝐴+𝑅≤𝜇. (2.4) It is clear that 𝑃 is an invariant set with respect to system (2.1) and any solution of system (2.1) initiated elsewhere will finally converge to the set.

3. Stability Analysis

3.1. The Basic Reproductive Number

In many epidemiological models, there exists a threshold condition, which determines whether an infection will be eliminated or become endemic. The basic reproduction number represents the average number of secondary infections caused by a single infected individual in susceptible population. For simplicity, denote 𝑎=𝜇+𝜇_𝐷+𝜐,𝑏=𝜇+𝜇^′_𝐷+𝑟₁,𝑐=𝜇+𝜇_𝐷^′′+𝑟₂,𝑎^′=𝑘𝑚𝜐+𝑘₁𝑝𝑎,𝑏^′=𝑘𝑚𝜐𝜎+𝑘₁1−𝑝𝑎(𝑏+𝜎)+𝑘₁𝑝𝑎𝜎+𝑘𝜐(1−𝑚)(𝑏+𝜎), (3.1) and further 𝑋=𝐴𝜇+𝜔𝛽₁𝑘𝑐(𝑏+𝜎)+𝛽₂𝑐𝑎^′+𝛽₃𝑏^′,𝑌=𝜇𝜔𝑟₁𝑐𝑎^′+𝑟₂𝑏^′,𝑍=𝜇𝜇+𝜔𝑎(𝑏+𝜎)𝑐𝜇+𝑘₁+𝑘. (3.2) Note that 𝑍−𝑌=𝜇𝜇+𝜔𝑎(𝑏+𝜎)𝑐𝜇+𝑘₁+𝑘−𝜇𝜔𝑟₁𝑐𝑎^′+𝑟₂𝑏^′>0. (3.3)

Computing the basic reproduction number yields 𝑅₀=𝑋+𝑌𝑍=𝐴𝜇+𝜔𝛽₁𝑘𝑐(𝑏+𝜎)+𝛽₂𝑐𝑎^′+𝛽₃𝑏^′+𝜇𝜔𝑟₁𝑐𝑎^′+𝑟₂𝑏^′𝜇𝜇+𝜔𝑎𝑐(𝑏+𝜎)𝜇+𝑘₁.+𝑘 (3.4) Obviously, 𝑅₀ is increasing with respect to 𝛽₁, 𝛽₂, 𝛽₃, and 𝜔. Note also that 𝛽₂, 𝛽₃ are less than 𝛽₁ after treatment. As a result, for sakes of deceasing 𝑅₀, it is important to increase the detected rate for disease control.

Evidently, system (2.1) has a disease-free equilibrium 𝑃₁(𝐴/𝜇,0,0,0,0,0). When 𝑅₀>1, it has yet a unique endemic equilibrium 𝑃₂(𝑆^∗,𝐸^∗,𝐼^∗_𝑈,𝐼^∗_𝐺,𝐼^∗_𝑅,𝑅^∗)∈𝑃⁰, the interior of 𝑃. Routine calculations show that 𝐼^∗_𝑈=𝜇𝐴𝜇+𝜔𝑘𝑐(𝑏+𝜎)(𝑋+𝑌−𝑍)𝑋()=()𝑅𝑍−𝑌𝜇𝑘𝑐𝐴𝜇+𝜔𝑏+𝜎₀−1(),𝑋𝑍𝑍−𝑌 (3.5) and 𝑆^∗, 𝐸^∗, 𝐼^∗_𝐺, 𝐼^∗_𝑅, and 𝑅^∗ could be obtained accordingly.

3.2. Proofs of Main Theorems

Theorem 3.1. The disease-free equilibrium 𝑃₁ of system (2.1) is globally asymptotically stable if 𝑅₀≤1 and is unstable if 𝑅₀>1.

Proof. In the region 𝑃, we take the Lyapunov function 𝐿𝑆,𝐸,𝐼_𝑈,𝐼_𝐺,𝐼_𝑅,𝑅=𝑆−𝑆₀−𝑆₀𝑆ln𝑆₀𝜔+𝐸+ℎ𝜔+𝜇𝑅+𝜇𝑐𝐼𝜇+𝜔^𝑅+ℎ₁𝐼_𝐺+ℎ₂𝐼_𝑈, (3.6) where 𝑆₀=𝐴𝜇,ℎ=𝐴𝛽₃𝜔+𝜇+𝜇𝜔𝑟₂,ℎ₁=𝐴𝛽₂𝜇+𝜔𝑐+𝜎ℎ+𝜇𝜔𝑐𝑟₁𝜇𝑐(),ℎ𝜇+𝜔𝑏+𝜎₂=𝐴𝛽₁𝜇+𝜔(𝑏+𝜎)𝑐+ℎ(𝑏+𝜎)𝑣(1−𝑚)𝜇𝑎()𝑐+𝜇+𝜔𝑏+𝜎𝜇𝜔𝑐𝑟₁𝑣𝑚+ℎ𝜎𝜐𝑚+𝐴𝛽₂𝜇+𝜔𝑐𝜐𝑚𝜇𝑐().𝜇+𝜔𝑏+𝜎 (3.7) It can be verified that 𝐿(𝑆,𝐸,𝐼_𝑈,𝐼_𝐺,𝐼_𝑅,𝑅)≥0 and 𝐿(𝑆,𝐸,𝐼_𝑈,𝐼_𝐺,𝐼_𝑅,𝑅)=0 if and only if (𝑆,𝐸,𝐼_𝑈,𝐼_𝐺,𝐼_𝑅,𝑅)=(𝐴/𝜇,0,0,0,0,0)=𝑃₁. This means that the function is positive definite.

Routine calculation shows that 𝑑𝐿=𝑆𝑑𝑡2−⁰𝑆−𝑆𝑆₀𝐴+𝜇+𝑘₁+𝑘𝑍(𝑋+𝑌−𝑍)𝐸. (3.8) When 𝑅₀≤1, we have 𝑑𝐿𝑆𝑑𝑡=𝐴2−⁰𝑆−𝑆𝑆₀+𝜇+𝑘₁+𝑘𝑍(𝑋+𝑌−𝑍)𝐸≤0. (3.9) It is clear that 𝑑𝐿/𝑑𝑡=0 if and only if either 𝑆=𝑆₀, 𝐸=0 when 𝑅₀<1 or 𝑆=𝑆₀ when 𝑅₀=1. The maximum invariant set in the region {(𝑆,𝐸,𝐼_𝑈,𝐼_𝐺,𝐼_𝑅,𝑅)∈𝑃,𝑑𝐿/𝑑𝑡=0} is {𝑃₁} when 𝑅₀≤1. It follows from the LaSalle invariant set theory that 𝑃₁(𝐴/𝜇,0,0,0,0,0) is globally asymptotically stable when 𝑅₀≤1.

The Jacobian matrix of system (2.1) at 𝑃₁(𝐴/𝜇,0,0,0,0,0) is computed to be 𝐽𝑃₁=𝐴₁𝐴₂𝐴₃𝐴₄, (3.10) where 𝐴₁=−𝜇0−𝛽₁𝐴𝜇0−𝜇+𝑘+𝑘₁𝛽₁𝐴𝜇,𝐴0𝑘−𝑎₂=−𝛽₂𝐴𝜇−𝛽₃𝐴𝜇0𝛽₂𝐴𝜇𝛽₃𝐴𝜇𝑤,𝐴000₃=0𝑘₁𝑝𝑚𝜐0𝑘₁𝜐(),𝐴1−𝑝1−𝑚000₄=−()𝑟𝑏+𝜎00𝜎−𝑐0₁𝑟₂−.𝜔+𝜇 (3.11) It can be obtained that |𝐽(𝑃₁)|=𝑍−(𝑋+𝑌). When 𝑅₀>1, we have |𝐽(𝑃₁)|<0. Then, it can be seen that 𝐽(𝑃₁) has at least an eigenvalue with positive real part, thus 𝑃₁ is unstable if 𝑅₀>1. By the Lyapunov method, equilibrium 𝑃₁ is unstable.

Remark 3.2. When 𝑅₀>1, we see that 𝐽(𝑃₁) has at least one eigenvalue with positive real part and at least one eigenvalue with negative real part, hence 𝑃₁ is a saddle point with a stable manifold and an unstable manifold. Solutions starting sufficiently close to 𝑃₁ must leave a neighborhood of 𝑃₁, except those on invariant 𝑇-axis which approach 𝑃₁ along the 𝑇-axis.

Theorem 3.3. System (2.1) has a forward bifurcation at 𝛽₁=𝛽^∗₁, and the endemic equilibrium 𝑃₂ is locally asymptotically stable if 0<𝑅₀−1≪1.

Proof. We proceed along the lines as Theorem  4.1 and Remark  1 in [6]. Suppose that 𝛽₁ is chosen to be a bifurcation parameter. By solving 𝑅₀=1 with variable 𝛽₁, we obtain 𝛽^∗₁=𝜇𝑎𝑐𝜇+𝜔(𝑏+𝜎)𝜇+𝑘₁+𝑘()−𝐴𝛽𝐴𝑘𝑐𝜇+𝜔𝑏+𝜎𝜇+𝜔₂𝑐𝑎^′+𝛽₃𝑏^′+𝜇𝜔𝑟₁𝑐𝑎^′+𝑟₂𝑏^′().𝐴𝑘𝑐𝜇+𝜔𝑏+𝜎 (3.12) When 𝛽₁=𝛽^∗₁, we have |𝐽(𝑃₁)|=0. When 𝛽₁>𝛽^∗₁, we have |𝐽(𝑃₁)|<0. When 𝛽₁<𝛽^∗₁, we have |𝐽(𝑃₁)|>0. Denote 𝐽(𝑃₁,𝛽^∗₁)=𝐽(𝑃₁)|_𝛽1=𝛽^∗₁.

First, we will show that 0 is a single eigenvalue of 𝐽(𝑃₁,𝛽^∗₁). To this end, let 𝐹(𝜆)=det(𝐽(𝑃₁,𝛽^∗₁)−𝜆𝐼) and 𝐽(𝑃₁,𝛽^∗₁)−𝜆𝐼=(𝑎_𝑖𝑗(𝜆)), 𝑖,𝑗=1,…,6. Denote the derivative of 𝐹(𝜆) by 𝐹^′(𝜆) and the algebraic complement of 𝑎_𝑖𝑖(𝜆) by 𝐴_𝑖𝑖(𝜆) for 𝑖,𝑗=1,…,6. Simple calculation gives 𝐹^′(0)=−𝐴₁₁(0)+𝐴₂₂(0)+𝐴₃₃(0)−𝐴₄₄(0)+𝐴₅₅(0)+𝐴₆₆(0)>0. (3.13) As a result, 0 is simple eigenvalue of 𝐽(𝑃₁,𝛽^∗₁). It can be seen that all the other eigenvalues of 𝐽(𝑃₁,𝛽^∗₁) have negative real parts.

Denote, respectively, the right and the left eigenvectors of 𝐽(𝑃₁) corresponding to the zero eigenvalue by 𝐿=𝐿₁,𝐿₂,𝐿₃,𝐿₄,𝐿₅,𝐿₆^𝑇,𝑀=𝑀₁,𝑀₂,𝑀₃,𝑀₄,𝑀₅,𝑀₆. (3.14) It follows from 𝐽𝑃₁,𝛽^∗₁𝐿=0,𝑀𝐽𝑃₁,𝛽^∗₁=0, (3.15) that 𝐿₁=𝑌−𝑍𝜇²()𝐿𝑘𝑐𝜔+𝜇𝑏+𝜎³,𝐿₂=𝑎𝑘𝐿₃,𝐿₃𝐿≠0,₄=𝑎^′𝑘()𝐿𝑏+𝜎³,𝐿₅=𝑏^′()𝐿𝑘𝑐𝑏+𝜎³,𝐿₆=𝑐𝑟₁𝑎^′+𝑟₂𝑏^′𝑘𝑐()𝐿𝜔+𝜇𝑏+𝜎³,𝑀₁=0,𝑀₂𝑀≠0,₃=𝐴𝛽^∗₁𝑐𝜇+𝜔(𝑏+𝜎)+𝑚𝜐ℎ₃()𝑀𝜇𝑎𝑐𝜇+𝜔𝑏+𝜎²+(1−𝑚)𝜐(𝑏+𝜎)𝐴𝛽₃𝜔+𝜇+𝜇𝜔𝑟₂()𝑀𝜇𝑎𝑐𝜇+𝜔𝑏+𝜎²,𝑀₄=ℎ₃()𝑀𝜇𝑐𝜇+𝜔𝑏+𝜎²,𝑀₅=𝐴𝛽₃𝜔+𝜇+𝜇𝜔𝑟₂𝑀𝜇𝑐𝜇+𝜔²,𝑀₆=𝜔𝑀𝜔+𝜇², (3.16) where ℎ₃=𝐴𝑐𝛽₂𝜔+𝜇+𝐴𝛽₃𝜎𝜔+𝜇+𝜇𝜔𝑟₂𝜎+𝜇𝜔𝑐𝑟₁. (3.17) Assume that 𝑀₂>0 and 𝐿₃>0 which are always possible up to a sign change of 𝐿 and 𝑀. The conditions in [6] are satisfied. Furthermore, routine calculation shows that the quantities as in [6] are 𝑀₂𝐿₁𝛽^′₁𝐿₃+𝛽₂𝐿₄+𝛽₃𝐿₅𝑀<0,₂𝐴𝜇𝛽^′₁𝐿₃+𝛽₂𝐿₄+𝛽₃𝐿₅>0, (3.18) respectively. By applying [6], we show that system (2.1) has a forward bifurcation at 𝛽₁=𝛽^∗₁. Therefore, when 𝛽₁>𝛽^∗₁ (i.e., 𝑅₀>1) with 𝛽₁ close to 𝛽^∗₁, 𝑃₂ is locally asymptotically stable.

4. Simulations

We use TB data provided by the Antituberculosis Research Institute of Guangdong Province (ARIG). We collect TB data from 1992 to 2008 and estimate the trend of TB cases in the near future in Guangdong province. The World Bank Loaned China Infectious and Endemic Disease Control Project has been carried out since 2000. According to the project plan, Guangdong province would implement the DOTS with the coverage rate 100%, the discovery rate 80%, and the treatment rate 90% in 2010. In this section, we carry out some simulations to describe the dynamic behavior of system (2.1), and predict the tendency of the TB disease. Setting one year as unit time, first we determine the values of the parameters using the data from 1992 to 2002, and then, we verify the model using the data from 2003 to 2008.

Assume that 𝛽₁=0.00000022, 𝛽₂=0.00000002, 𝛽₃=0.00000006, which are adjusted from the literature [2, 5]. Furthermore, assume that 𝑘+𝑘₁=0.0052 based on the estimation that 10% of the latent cases will develop into infectious person in 20 years [9]. The detected rate in time 𝑘₁ and the treatment failure rate 𝜎 that we estimate are, respectively, 3.2×10⁻³ and 0.04. Among the detected infectious patients, about 8.32% are MDR-TB patients from our data, so we assume that 𝑚=𝑝=0.916. The other parameters in Table 1 can be estimated or derived from our data using the weighted average method. In Table 1, GHSY means Guangdong health statistical yearbook.

Table 1 

Model parameters for simulations.

Using the investigated data in 2002 in Guangdong province as the initial population sizes 𝑆(02),𝐸(02),𝐼_𝑈(02),𝐼_𝐺(02),𝐼_𝑅()(=[],02,𝑅02)41241970,33600080,158990,105250,5290,46810 (4.1) we predict the TB situation from 2003 to 2008 in Figure 2(a). ⋆ and ▲ represent predicted infection TB cases and actual TB cases, respectively. That the prediction basically coincides with the actual infection number demonstrates the ability of our model in describing the TB disease.

(a)

(b)

Figure 2 

Trend simulation and prediction of the infection cases in Guangdong, China. (a) Contrast diagram of predicted TB infection cases and actual TB infection cases. (b) Prediction result of TB infection cases.

To obtain simulation results in a more accurate manner, we adjust some parameters by using data from 1992 to 2008. Let 𝐴=950000, 𝜇=0.013, 𝜇^′_𝐷=0.0023, and 𝜇_𝐷^′′=0.02 in the following simulations. Using the other parameter values given in Table 1 and the data of 2008 𝑆(08),𝐸(08),𝐼_𝑈(08),𝐼_𝐺(08),𝐼_𝑅()(=[],08,𝑅08)40241000,38600000,98990,155200,25200,96800 (4.2) as the initiate state, we obtain the predicted data shown in Figures 2(b) and 3. We observe that the population of the undetected infected individuals slowly tends to a stable endemic equilibrium. The total TB infected number will be 3.28×10⁵ in 2015. This indicates a very severe situation, which concurs with the warning from the leading experts [10].

Figure 3 

Simulation of system (2.1) in Guangdong province, China. Disease persistence in the population. Using the parameters in Table 1, we obtain 𝑅0=3.5846>1.

Figure 4(a) is a graphical representation of total TB infection number under three different delay discovery rates. It is clear that only when 𝜐≥0.9, the infection number will keep relatively stable.

(a)

(b)

(c)

(d)

Figure 4 

Simulations with different parameters. (a) Different discovery rates. (b) Different detected rates in time. (c) Different treatment rates of the MDR-TB patients. (d) All TB infection cases are detected in time.

Figure 4(b) is a graphical representation of total TB infection number with different 𝑘₁’s. It shows the influence of the detected rate in time. If the delay rate 𝑘/(𝑘+𝑘₁) is about 80%, the total infection number rapidly raises. On the contrary, if the detected rate in time 𝑘₁=0.0042, the amount of the infection number will basically remains stable from 2014. It can be seen that the detected rate in time is an important factor for disease control.

Figure 4(c) is a graphical representation of total TB infection number with different treatment rate of the MDR-TB patients. If 𝑟₂=0.8, the total infection number will be smaller than 3×10⁵ in 2012.

5. Discussions and Conclusions

In this work, combining stability analysis with bifurcation identification, we have investigated the dynamical behavior of a TB model and taken into account two features, the undetected TB cases and MDR-TB cases. The model had an influence up on TB prevention and control. We found that there are two possible steady states: the healthy and the infected equilibria. We proved the global asymptotic stability of the disease-free equilibrium when the basic reproduction number is equal to or less than unity. The local asymptotic stability of the endemic equilibrium has been also established when the basic reproduction number is bigger than unity. A more delicate analysis based on the center manifold theory indicated that a forward bifurcation existed under a proper condition. Forward bifurcation showed that there existed a stable endemic equilibrium when the disease-free equilibrium become unstable. Furthermore, this work does not give proofs of global asymptotical stability of the endemic equilibrium.

The number of the MDR-TB cases in China is about one third as that of the world. As a province with the one-tenth TB cases in China, Guangdong province suffers much from the TB disease. From Figure 3, we can find the undetected infected patients slowly increase. A primary reason for this is that the detection ratio has been improved thanks to the implementation of the World Bank Loaned China Infectious and Endemic Disease Control Project. However, due to the facts that the total population of residence is huge in Guangdong (about 100 millions) and more than 20 percent of residents are temporary immigrants from the rural inland, the number of infected people will slowly go up. The huge population of rural migrants makes the in-time detection more difficult. We have analyzed the data collected by ARIG and carried out simulations based on our model. First, the parameters of the model were identified, which justifies the important role that undetected cases play in the model. Second, we predicted the trend of number of TB infection cases in 20 years and found that it was difficult to eliminate the TB disease under current conditions. Third, we analyzed the possibility of controlling the disease by improving the detected rates and other factors. In particular, when every patient can be detected in time, that is, 𝑘=0, 𝑘₁=0.0052, and the detected patients have no infectiousness, then we can obtain 𝑅₀=0.1879<1 even if 𝜐=0.6 and 𝑟₂=0.6. Under these assumptions, the solutions of the system will tend to the disease-free equilibrium 𝑃₁, and the TB disease will finally die out as indicated in Figure 4(d). This once again shows the importance of the detection rate in TB control.

Acknowledgments

The authors are grateful to the referees for the valuable suggestions and numerous comments which improve the quality of the paper. This work was supported by China’s National Science and Technology Special Projects during the 11th Five-Year Plan Grant no. 2008ZX10003-007 and the National Natural Science Foundation of China under Grants nos. 60925013 and 30973981.