key: cord-0975643-sy4gwigp
authors: Cai, Yongli; Zhao, Shi; Niu, Yun; Peng, Zhihang; Wang, Kai; He, Daihai; Wang, Weiming
title: Modelling the effects of the contaminated environments on tuberculosis in Jiangsu, China
date: 2020-09-16
journal: J Theor Biol
DOI: 10.1016/j.jtbi.2020.110453
sha: bc2223ef131ab11581e93bc506c7a4ea46867207
doc_id: 975643
cord_uid: sy4gwigp

Tuberculosis (TB) is still an important public health issue in Jiangsu province, China. In this study, based on the TB transmission routes and the statistical data of TB cases in Jiangsu, we formulate a novel TB epidemic model accounting for the effects of the contaminated environments on the TB transmission dynamics. The value of this study lies in two aspects: Mathematically, we define the basic reproduction number, [Formula: see text] , and prove that [Formula: see text] can be used to govern the threshold dynamics of the model. Epidemiologically, we find that the annual average [Formula: see text] is [Formula: see text] and TB in Jiangsu is an endemic disease and will persist for a long time. As a consequence, our work suggests that, in order to control the TB in Jiangsu efficiently, we must decrease the virus shedding rate or increase the recovery rates, and increase the environmental clearance rate.

1 Tuberculosis (TB) is an ancient and chronic infectious disease, which is caused by infection 2 with the Mycobacterium tuberculosis (MTB) [1] . Normally, the TB bacteria are put into the air 3 when a person with TB disease, and TB patients are mainly transmitted by droplets produced by potentially infective respiratory droplets. Transmission may also occur through fomites in the patients [10] , include alcoholism [9] and diabetes mellitus (three-fold increase) [11] . 23 Currently, there are approximately 95% of the estimated 8 milion new cases of TB occuring 24 in developing countries each year, and two-thirds of which appears in India (27%), China (9%), 25 Indonesia (8%), Philippines (6%), Pakistan (5%), Nigeria (4%), Bangladesh (4%) and South Africa 26 (3%), where 80% occur among people between the ages of 15 to 59 years. And only 6% of global 27 cases occur in Europe (3%) and the Americas (3%) [2] . 28 In recent years, the Chinese government has increased its investment in public health, and 29 the laws and regulations for disease prevention and control has been constantly improved, which 30 provides an important basic guarantee for coping with major public health emergencies, preventing 31 infectious diseases and ensuring the health of the people, thus effectively controlling major infectious 32 diseases. In particular, the incidence of tuberculosis has decreased significantly. According to the . 37 It was known that 80% of TB exists in rural areas, particularly in north and north-western 38 regions with low socioeconomic status in China [13] . But in the past 20 years in Jiangsu province, 39 China (see Fig. 1 for the location of Jiangsu province in China), one of the most developed areas in 40 China in economy, technology and culture and the total output is one of the largest in the nation, 41 TB ranked first in the number of notifiable B infectious diseases [14] . In 2011, Jiangsu had 39,589 42 existing TB, and in 2017, 28,402 existing TB and in 2018, there is 26,506 incident TB. In particular, 43 the incidence of tuberculosis has decreased significantly (see Fig. 3 (a) for more details). However, 44 the situation of prevention and control of TB in Jiangsu province is still very serious. 45 It is worthy to notice that mathematical models have played a key role in the formulation of TB 46 control strategies. Waaler et al. [15] introduced the first mathematical model for TB in ordinary 47 differential equations. The simplest TB transmission models include classes of susceptible, exposed, 48 and infectious individuals, and hence, they are known as the SEI models. Of course, there are more The main focus of this paper is to investigate how contaminated environments affect TB dy-68 namics through studying the threshold dynamics of a general TB model. And the rest of the paper 69 is organized as follows. In Section 2, we formulate the model in details. In Section 3, we give the 70 dynamics analysis of the model, we introduce the basic reproduction number R 0 and prove that 71 R 0 can be used to govern the threshold dynamics of the model. In Sections 4, we give the TB 72 epidemics in Jiangsu, China via numerical simulations. In the last section, Section 5, we provide a 73 brief discussion and the summary of the main results. Suppose that the total population individuals N (t) divide into susceptible S(t), infectious but 76 not yet symptomatic, i.e., pre-symptomatic T e (t), infectious with symptoms T (t), and recovered 77 R(t). We further consider the TB virus concentration in environment as W (t), which is the density 78 of pathogen of the contaminated environments including door handles, towels, handkerchiefs, toys, 79 utensils, bed and toilet seat, bathroom washbasin tap lever, bathroom ceiling-exhaust louvre etc. 80 at time t. And our model involves two typical transmissions: one is the direct transmission between 81 susceptible S(t) and infected individuals (including pre-symptomatic T e (t) and symptomatic T (t)) 82 with rate of β 1 (t); the other is the indirect transmission to susceptible individuals and infected indi-83 viduals (i.e., T e (t) and T (t)) by contaminated environments W (t) with rate of β 2 (t). A seasonality 84 in the long-term patterns of TB incidences time series can be observed evidently (see Fig. 3 ), and 85 there is a growing awareness that seasonality can cause population fluctuations ranging from annual 86 cycles to multiyear oscillations [21] , and hence we assume that the transmission rates β 1 and β 2 (t) Thus we can establish the following TB epidemic model involving five ordinary differential 90 equations: and N = S + T e + T + R.

The meanings of each variables and parameters in model (2.1) are as follows.

• µ: the per capita natural mortality rate; • υ: reactivation rate the pre-symptomatic infectious;

• γ 1 : the recovery rate of the pre-symptomatic infectious;

• γ 2 : the recovery rate of the symptomatic infectious;

• ρ: the rate from recovered to susceptible.

• δ: disease-related death;

• α: the virus shedding rate from symptomatic infected individuals;

• c: the clearance rate of the virus in the environments.

For the sake of conveniently analysis, set u(t) = (u 1 , u 2 , u 3 , u 4 , u 5 ) = (T e , T, W, R, S). Then 109 model (2.1) with (2.2) is equivalent to the following system:

with initial conditions u(0) = u 0 , u i (0) = u 0 i (i = 1, 2, · · · 5).

Theorem 2.1. Model (2.3) has a unique and bounded solution with the initial value u 0 ∈ R 5 + , i.e,

Proof. Considering the non-negativity of I, i.e., I ≥ 0, the total population N (t) can be determined 111 by the following model:

(2.4)

It is easy to see that the linear differential equation dN dt = Λ−µN has a unique equilibrium N * = Λ µ , 113 which is globally asymptotically stable. The comparison principle implies that lim t→∞ N (t) ≤ Λ µ .

Then ∀ ε > 0, there exists a t 0 > 0 such that 115 N (t) ≤ N * + ε, t > t 0 .

(2.5)

It follows from the comparison principle that lim t→∞ u 3 (t) ≤ αΛ cµ . This completes the proof.

Let Φ A(·) (t) be the fundamental solution matrix of equation

continuous, cooperative, irreducible and ω-periodic n × n functional matrix. Let r(Φ A(·) (ω)) be the 118 spectral radius of Φ A(·) (ω). 

Obviously model (2.3) admits a disease free equilibrium (DFE) E 0 = 0, 0, 0, 0, Λ µ . Then

Let Y (t, s) (t ≥ s) be the evolution operator of the linear ω-periodic system

Following the method established by Wang and Zhao [28] , let φ(s) be ω-periodic in s and the initial distribution of infectious individuals. So F (s)φ(s) is the rate of new infections produced by the infected individuals who are introduced at time s. When t ≥ s, Y (t, s)F (s)φ(s) gives the distribution of those infected individuals who are newly infected by φ(s) and remain in the infected compartments at time t. Naturally,

is the distribution of accumulative new infections at time t produced by all those infected individuals 128 φ(s) introduced at time previous to t. 129 Let C ω be the ordered Banach space of all ω-periodic functions from R to R 3 , which is equipped with the maximum norm · and the positive cone

Then we can define a linear operator L implies that

which is called the next infection operator, and the spectral radius of L is defined as the basic 130 reproduction number:

(3.2)

In order to characterise R 0 , we introduce the linear ω-periodic system

. Following the general calculation procedure in [28, Theorem 2.1], the basic reproduction number 133 R 0 is the unique solution of r(W (ω, 0, λ)) = 1.

Following Wang and Zhao [28], we can obtain the relation between R 0 and r(W (ω, 0, λ)) = 1 135 shown in the following lemma.

Lemma 3.1. [28, Theorem 2.2] The following statements are valid:

Remark 3.2. In the special case of β 1 (t) = β 1 and β 2 (t) = β 2 , the basic reproduction number R 0 of

Here Proof. Consider an auxiliary system

It follows Lemma 2.2 that there exits a positive ω-periodic function

Choose t 1 > t 0 and a small number

where T is the transposition of the vector.

It follows from R 0 < 1 that r(Φ F −V (ω)) < 1. Since r(Φ F −V +εM (ω)) is continuous for all small 161 ε, we can choose ε > 0 small enough such that r(Φ F −V +εM (ω)) < 1. Hence, we get p < 0. It follows In the following, we attempt to explore the uniform persistence of model (2.1) when R 0 > 1. Define

Lemma 3.5. X and X 0 are positively invariant.

Proof. From the fifth equation of model (2.3), we can derive that

By [29, Theorem 4.1.1] as generalized to nonautonomous systems, the irreducibility of the cooperative matrix

Thus, X and X 0 are positively invariant.

Clearly, ∂ X 0 is relatively closed in X.

Let P : X → X be the Poincaré map associated with model (2.3), that is,

From Theorem 2.1, we know that P is a dissipative point on R 5 + . Thus P admits a global attractor, 170 which attracts every bounded set in R 5 + . We then introduce the following lemma.

Lemma 3.6.

where d(P m (u 0 ), E 0 ) represents distance between P m (u 0 ) and E 0 .

Proof. Since R 0 > 1, Lemma 3.1 implies that r(Φ F −V (ω)) > 1. Then we can choose ǫ > 0 small

(3.7)

By the continuity of the solutions with respect to the initial values, for ǫ > 0, there exists a σ * > 0 such that ∀ u 0 ∈ X 0 with u 0 − E 0 ≤ σ * , we can obtain

We proceed by contradiction to prove that

If not, we can get lim sup

Without loss of generality, we can assume that d(P m (u 0 ), E 0 ) < σ * for all m ≥ 0. Then we can get

For any t ≥ 0, let t = mω + t ′ , where t ′ ∈ [0, ω] and m = t ω , which is the greatest integer less than or equal to t ω . Then,

176 Consider the following auxiliary system

.

It follows from Lemma 2.2 that there exists a positive ω-periodic functionv(t)

By the comparison principle, we have

Then (u 1 (t), u 2 (t), u(t)) T → ∞ as t → ∞, a contradiction. This completes the proof.

Define

Lemma 3.7. P is uniformly persistent with respect to (X 0 , ∂ X 0 ).

Proof. Now we first prove 181 M ∂ = {(0, 0, 0, 0, u 5 ) ∈ X, u 5 ≥ 0}.

(3.9)

Noting that

we only need to prove that M ∂ ⊆ {(0, 0, 0, 0, u 5 ) ∈ X, u 5 ≥ 0}.

It suffices to prove that for any u 0 ∈ M ∂ , we have u i (mω) = 0, ∀ m ≥ 0, i = 1, 2, 3, 4. If it is not true, there exists an m 1 ≥ 0, such that (u 1 (m 1 ω), u 2 (m 1 ω), u 3 (m 1 ω), u 4 (m 1 ω)) T > 0. Thus (3.6) implies that 

Thus, if u 0 ∈ {(0, 0, 0, 0, u 5 ) ∈ X, u 5 ≥ 0}, then u 0 ∈ M ∂ , which contradicts with u 0 ∈ M ∂ . Hence,

M ∂ ⊆ {(0, 0, 0, 0, u 5 ) ∈ X, u 5 ≥ 0}, which implies that (3.9) holds. Clearly, E 0 is the only fixed with u * 5 (0) = u * 5 (nω) = 0, n = 1, 2, 3 · · · , where a * (t) = µ + β 1 (t) u * 1 N + β 2 (t)u * 3 . Therefore, we have 

where a ∈ [0, 1].

To find endemic equilibrium, we make the substitution x = u 2 N . Then,

(3.12)

Since Λ − µN − δu 2 = 0, By using some algebraic computations, we can obtain

, where x is a positive real root of the following equation:

If R 0 ≤ 1, then R dir 0 ≤ 1, C ≥ 0 and B > 0. It follows that Eqn. (3.13) has no positive real root. If R 0 > 1, equation (3.13) has a unique real root x * :

Hence model (2.3) has a unique endemic equilibrium E * = (u * 1 , u * 2 , u * 3 , u * 4 , u * 5 ) with

(3.14)

In a special case of δ = 0, if R 0 > 1, model (2.3) has a unique endemic equilibrium Fig. 3(a) , we show the local TB epidemic from 2009 to 2018. 207 We can observe that a substantially decreasing trend in the TB incidences, dropped from some 4000 208 cases per month in 2009 to some 2000 cases per month in 2018.

The relative infectivity can be (preliminary) quantified by using the approach in Fine and [33, 34, 35] . The relative infectivity of TB, Q t , is quantified in Fig. 3(b) .

In Fig. 3(c) , we show the strong seasonality [31] in the TB incidence time series in the annualised 216 epidemic curves. And in Fig. 3(d) , we show the annualised relative infectivity to show the seasonality 217 in the TB infectivity across years. We can find that the (relative) infectivity in February appeared 218 to be dominant (or the highest) across different months. . 248 We attempt to find the MLEs of both R dir 0 and R ind 0 (t) within biologically and clinically rea- 

First of all, we show the parameters' values used for the numerical simulation and sensitivity 255 analysis for model (2.1) in Table 1 . to be strictly less than one. We can find that the estimated annual trends in R ind 0 (t) (Fig 4(b) ) are 259 consistent with the patterns in the preliminary descriptive "infectivity" in Fig 3(d) . We estimated 260 that the R ind 0 in February is larger than one (estimated to be 6.9 in Fig 5(b) ), whereas those in 261 other months are likely to below one.

In Fig. 4(c) , we give the fitting results and projection to 2019. We find that a decreasing trend 263 in the fitting results, which matched the trends in the observed incidence data. The projection is 264 also likely to maintain these trends in 2019.

It is should be noted that the term β 2 is considered as a time-varying parameter such that the 266 indirect transmission could also be time-varying. Hence, the R ind 0 (t) reconstructed in Fig. 4(b) is 

Since the numerical simulation is conducted with the same complexity in the model structure, 278 we can directly study the goodness-of-fit (in term of the likelihood) and the fitting errors. The in a Normally distribution (Fig 6(a)-(b) ). And the mean percentage error is close to zero (Fig 6(c) ). with changes in the epidemiological parameters. Fig. 8(a) and Fig. 8(b) show the trend changes in 

The main focus of this study is to investigate the effects of the contaminated environments on 297 the TB transmission dynamics in Jiangsu, China analytically and numerically. Mathematically, 298 we define the basic reproduction number R 0 (cf. (3.2 ), and prove that R 0 can be used to govern 299 the threshold dynamics of the model: if R 0 < 1, the unique DFE is globally asymptotic stable 300 (cf. Theorem 3.4); while R 0 > 1, there is at least one positive periodic solution and TB will 301 persist uniformly (Theorem 3.9). Epidemiologically, we show that the cost of the contaminated 302 environments affect the transmission dynamics of TB in Jiangsu, China in the following aspects: 303 (i) Based on the monthly TB incident cases counted by the Jiangsu CDC (cf. Fig. 3(a) ), the 304 TB incidence time series has strong seasonality (cf. Fig. 3(c) ). And the annual average 305 R 0 is 1.13 > 1, then from Theorem 3.9, we can conclude that the TB in Jiangsu persists 306 under current circumstances. That is, the TB becomes an endemic disease and will persist in 307 Jiangsu for a long time. And there is a long way to go to achieve the world-wide goal towards 308 elimination of tuberculosis by 2050. (ii) The annualised relative infectivity (cf. Fig. 3(d) ) shows that the relative infectivity in Febru-310 ary is dominant across different months. This is consistent with the estimation of that R ind 0 311 in February is 6.9 (cf. Fig 5(b) ), whereas those in other months are likely to below one. This Table 1 having a random perturbation with a coefficient of variation of 0.2. Table 1 . less than one, we remark that the TB epidemics are likely to be controlled, in term of R 0 < 338 1, by effectively controlling the indirect transmission path. Also, Fig. 8 indicates that the 339 control measures reduce the R ind 0 (t) or R dir 0 could effectively decrease the number of TB cases.

This could be achieved by providing timely and effective treatment and maintaining a low 341 environment contamination level. We further find that a lower virus shedding rate, α, and a 342 higher environmental clearance rate, c, lead to low level of environment contamination (Fig. 8 ) 343 and R 0 (Fig. 7) , which could control the TB epidemic efficiently.

(iv) From Fig 7, the TB transmissibility and number of cases are positively associated with the 345 effective transmission rates β 1 and β 2 , as well as the virus shedding rate α. The effective 346 control efforts are suggested to focus on reducing the β 1 , β 2 and α. We also find that the 347 recovery rates, γ 1 and γ 2 , and environmental clearance rate, c, are negatively associated with 348 the TB transmissibility and number of cases. Thus, increasing γ 1 and γ 2 are also likely to 349 control the TB epidemics. Obviously, increasing close-contact distance and wearing protective 350 mask can efficiently increase γ 1 , γ 2 and therefore useful to control the spread of TB in Jiangsu,

China.

It is worthy to note that, more and more evidences proved that the respiratory infections are 

Control strategies for tuberculosis epidemics: 378 new models for old problems

Air, surface environmental, and personal pro-381 tective equipment contamination by severe acute respiratory syndrome coronavirus

CoV-2) from a symptomatic patient

Criteria for the control of drug-resistant tuberculosis

Modelling the effects of contaminated environments on 386 HFMD infections in mainland China

Modelling seasonal HFMD infections with the effects of 388 contaminated environments in mainland China

Toilets dominate environmental detection of SARS-CoV-2 390 virus in a hospital. medRxiv

World Health Organization. Modes of transmission of virus causing COVID-19: impli-392 cations for IPC precaution recommendations

Tuberculosis: Disease of the past, disease of the present

Convergence of the tuberculosis and diabetes epidemics: renewal of old ac-398 quaintances

Survey of the epidemic situation of notifiable infectious diseases in China

The reported tuberculosis cases in Jiangsu province

The use of mathematical models in the study of the 405 epidemiology of tuberculosis

Prospects for worldwide tuberculosis control 407 under the WHO DOTS strategy

Potential public health impact of new tuberculosis vaccines

Quantifying the intrinsic transmission dynamics of tuberculosis

A model for tuberculosis with exogenous 413 reinfection. Theor

Rate of reinfection tuberculosis after suc-415 cessful treatment is higher than rate of new tuberculosis

A tuberculosis model with seasonality

The modelling of hand, foot, and mouth disease in contami-420 nated environments in

A periodic epidemic model in a patchy environment

On the definition and the computation of 427 the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations

Reproduction numbers and sub-threshold endemic 430 equilibria for compartmental models of disease transmission

Mathematical epidemiology of infectious diseases

Threshold dynamics for compartmental epidemic models in peri-434 odic environments

Monotone dynamical systems: an introduction to the theory of competitive and 436 cooperative systems

Dynamical systems in population biology

Measles in england and wales: an analysis of factors underlying 440 seasonal patterns

Modelling the large-scale yellow fever outbreak in luanda

angola, and the impact of vaccination

Different epidemic curves for severe acute respiratory syndrome 444 reveal similar impacts of control measures

Estimating individual and household reproduction numbers in an emerging epidemic

Simple framework for real-time forecast in 448 a data-limited situation: the Zika virus (ZIKV) outbreaks in brazil from 2015 to 2016 as an 449 example

Statistical inference for partially observed markov 451 processes via the R package pomp

Modelling the skip-and-resurgence of Japanese encephali-453 tis epidemics in Hong Kong

Plug-and-play inference for disease dynamics: measles in large 455 and small populations as a case study

Modeling the spread of middle east respiratory syndrome 457 coronavirus in Saudi Arabia

Inference for nonlinear dynamical systems

Monte carlo profile confidence 461 intervals for dynamic systems

R: A language and environment for statistical computing

A mathematical model to study the 2014-464 2015 large-scale dengue epidemics in Kaohsiung and Tainan cities in Taiwan

Prevention and control of Zika as a mosquito-borne and 467 sexually transmitted disease: A mathematical modeling analysis

Implication of vaccination against dengue for Zika outbreak

Modelling weekly vector control against dengue in the 471

Guangdong Province of China

Community transmission of severe acute respiratory 473 syndrome coronavirus 2

Early transmission dynamics in wuhan, china, of novel 475 coronavirus-infected pneumonia

Clinical features of patients infected with 2019 novel 477 coronavirus in Wuhan

Multi-route respiratory infection: when a transmission route 479 may dominate. medRxiv

Global burden of drug-resistant tuberculosis in 481 children: a mathematical modelling study

Global burden of latent multidrug-resistant 483 tuberculosis: trends and estimates based on mathematical modelling

Ambient air pollution exposures and risk of drug-resistant 486 tuberculosis

Long-term exposure to ambient air pollution 488 and mortality in a Chinese tuberculosis cohort

Ambient particulate air pollution and daily mortality in 490 652 cities

The authors would like to thank the editor and the referees for their helpful comments. Y. Cai