=Paper=
{{Paper
|id=Vol-1904/paper22
|storemode=property
|title=Study of a singularly perturbed tuberculosis model
|pdfUrl=https://ceur-ws.org/Vol-1904/paper22.pdf
|volume=Vol-1904
|authors=Elena Tropkina,Elena Shchepakina
}}
==Study of a singularly perturbed tuberculosis model ==
https://ceur-ws.org/Vol-1904/paper22.pdf
Study of a singularly perturbed tuberculosis model
E. Tropkina1, E. Shchepakina1
1
Samara National Research University, 34 Moskovskoe Shosse, 443086, Samara, Russia
Abstract
A detailed analysis of the dynamic model of the tuberculosis epidemic was conducted. It is shown that the dynamic model contains several
time scales and can be represented in a singularly perturbed form. With help of the integral manifolds theory and the reduction principle, the
reduction of the modeling system was justified and carried out. This approach allow us to replace the original system by another system of a
lower order on an integral manifold whose dimension is equal to that of the slow subsystem and which, at the same time, preserves the
essential properties of the original system. The conditions for the stabilization of the epidemiology based on the selection of necessary
treatment and preventive measures are determined.
Keywords: singular perturbations; integral manifold; order reduction; stability; epidemiology; tuberculosis
1. Introduction
It is well known that tuberculosis is a deadly disease, the fight against which is still relevant today. In 1882, Robert Koch,
along with the discovery of the tubercle bacillus, found that this disease is transmitted by aerogenic means [1]. Consequently,
people who frequently contact individuals with an active form of tuberculosis (the infectious stage of the disease) have a much
higher risk of infection. Majority of infected people remain latent carriers throughout their life. The average duration of the
latent period (the period of latent infection) ranges from several months to dozens of years. However, the risk of progression to
the active form of tuberculosis increases dramatically in the presence of concomitant infections that weaken the immune system.
In the absence of treatment for tuberculosis of respiratory organs, mortality is about 50%.
Currently, about 3 million people die from tuberculosis every year worldwide [2]. But in most cases, tuberculosis is curable.
The current methods of treating this disease require long-term courses of treatment (from six months to several years), the
violation of which often leads to the return of the disease and the development of drug resistance.
Treatment of drug-resistant tuberculosis is much more dangerous for the patient, less successful and more costly than
treatment of the disease caused by common strains of the pathogen. All this shows the necessity of developing a new approach
of identifying and treating patients with tuberculosis, as well as timely and effective treatment.
One of the most effective methods for solving such problems is the construction and study of a mathematical model
describing the processes of infection spread in the population, the development of the disease and the impact of anti-tuberculosis
measures. Based on this study we can determine effective measures to combat this dangerous phenomenon.
In the present work, a detailed analysis of the dynamic model of the tuberculosis epidemic based on the cluster approach has
been carried out. The presence of several time scales made it possible to apply the geometric theory of singular perturbations for
its qualitative study. This approach made it possible to determine the conditions for the stabilization of epidemiological
conditions based on the selection of the necessary treatment and preventive measures.
2. Cluster model of the tuberculosis epidemic
Mycobacterium drops of tuberculosis get into the air when coughing or sneezing infected people. Tuberculosis bacillus,
spreading by such drops, lives in the air for a short period of time (about two hours) and, therefore, it is believed that casual
contact with persons with active tuberculosis (infected individuals) rarely lead to the spread of the disease, and that most relapses
are the result of prolonged and close contacts with primary carriers of infection. Latently infected individuals become infectious
after some, usually long, time period. This period of transition to an active form of infection is called latent. Latent periods vary
from several months to dozens of years. Most infected individuals never go to the active form of tuberculosis. On the other hand,
the average length of the infection period is relatively short (several months). This indicator is decreasing in countries with
affordable treatment.
A common scheme for analyzing the spread of the tuberculosis epidemic is based on the division of the population into
certain classes. Consider the basic mathematical model of the spread of tuberculosis [2]:
dS I
dt S 1 S N ,
dE S I ( k r ) E E I R I ,
dt 1
N
1 2
N
3
N
(1)
dI I
k E E ( d r ) I ,
dt 2
N
2
dR I
r2 I r1 E R 3 R .
dt N
3rd International conference “Information Technology and Nanotechnology 2017” 117
Mathematical Modeling / E. Tropkina, E. Shchepakina
Model (1) was derived under assumption that the total population is divided into four epidemiological classes: susceptible (S),
i.e., uninfected, but susceptible to infection individuals, those in whose body the pathogens of tuberculosis have not yet
penetrated; carriers of latent infection (E), i.e., individuals in whose body tuberculosis pathogens are present, in equilibrium with
the immune system, such individuals are characterized by the absence of any external manifestations of the disease; infected
individuals (I) with clinical manifestations of tuberculosis caused by sufficiently extensive tissue damage as a result of the
activity of mycobacteria in their organisms; recovery individuals (R) who have completed treatment and recovered from the
disease. The parameter reflects the influx of young people into the model population The parameters 1 , 2 , 3 are the
transmission rates of tuberculosis infection for the respective classes; is natural mortality rate; k is the tuberculosis
progression rate; d is the tuberculosis induced mortality rate; r1 , r2 denote the treatment rates for latent class and infectious
class, respectively; N S E I R is the total population size. The mixing of classes in this model is homogeneous, that is,
there are no assumed differences between individuals while tuberculosis transmission depends on the rate of infection.
The basic reproduction number, one of the most important characteristics in mathematical biology, defined as the average
number of secondary cases produced by typical infected individuals, mainly in a susceptible population, is described by
1 k k
0HM Q0 ,
( d r2 ) ( k r1 ) (k r1 )
1
Q0 , d r2 ,
where Q0 is the number of secondary latent infected individuals produced by a typical infected individual during the mean
infectious period 1 , f k / (k r1 ) shows the probability of survival during the transition from the latent to the active
infectious stage.
It is assumed that only individuals who have frequent and prolonged interactions with infected people have a high risk of
tuberculosis infection. Newly infected individuals activate clusters (groups of individuals who come into regular and close
contact with people in an active form of tuberculosis), increasing the risk of tuberculosis infection for susceptible individuals of
each cluster [2].
Let us make a number of variables changes. We set a constant n be the average size of the cluster; the risk of infection with
tuberculosis in the cluster will be determined by the parameter . Further, the population of uninfected people in the cluster
will be given by N1 (t ) n I (t ), where N1 includes two subpopulations, namely, susceptible S1 and latent infected E1 , i.e.,
N1 (t ) n I (t ) S1 (t ) E1 (t ).
The population of persons who do not belong to the cluster at the time t is denoted as N2 . This population consists only of
.
the sensitive S2 and the latent infected E 2 individuals. The subpopulations of the sick are not included in this model for the
sake of simplicity. We assume that n k E2 (S2 / N2 ) individuals go to the class S1 per unit of time, while the individuals
n k E2 ( E2 / N2 ) go to the class E1 per unit time. In addition, since infected individuals are cured or die (with a speed I ),
that I is the rate at which clusters become inactive (or die). We suppose that n I S1 / N individuals returns to the class
S2 and n I E1 / N returns to the class E2 for the unit of time, respectively. Assuming a low level of distribution of
individuals with an active form of tuberculosis, we consider the case N1 N 2 , hence we can neglect the birth rate in the
population.
The above assumptions lead (1) to the following basic cluster model [2]:
dS1 S2
dt ( ) S1 N n k E2 ,
2
dE1 E2
S1 E1 n k E2 ,
dt N2
dI
k E2 I , (2)
dt
dS 2 S2
dt S 2 S1 N n k E2 ,
2
dE2 E
E1 ( k ) E2 2 n k E2 .
dt N 2
The basic reproduction number for system (2) is:
n k
c0 Q0 f .
( ) ( k )
Hence, Q0 n / ( ). is the expected number of infections produced by one infected individual in his cluster. Only a
part f k / ( k ) among the infected individuals will survive during the latent period.
3. Dimensionless model
Disease and dynamics of populations have characteristic time scales. Dynamics of influenza in humans is "super fast" at the
individual and community levels in comparison with the dynamics of the carriers of infection (people). This is because the
3rd International conference “Information Technology and Nanotechnology 2017” 118
Mathematical Modeling / E. Tropkina, E. Shchepakina
average life expectancy of an infected person is approximately 4000-8000 times the average duration of an influenza infection
and, therefore, the largest outbreaks of influenza occur in local communities before any significant demographic change can
occur (several months). Consequently, when studying the dynamics of the flu epidemic, two typical time scales are often used:
the time scale of the disease and the life span of the carrier of the infection [3].
Tuberculosis is usually described as a slow disease due to its long and varied distribution of the latent period and due to the
short and relatively narrow distribution of its infectious period [1]. Most latently infected with tuberculosis do not become
actively infected, i.e., there is no transition of latent form to active. Some become actively infected during a five-year period,
while others become active only after a longer period of time (perhaps decades). On the other hand, infected individuals remain
so for relatively short periods of time, in part because of the use of antibiotics (an average of six months). Since secondary
infections are formed from infected individuals, individuals with an active form of tuberculosis have a relatively short period for
possible infection of other people. Consequently, the infection of tuberculosis-susceptible individuals occurs on the same time
scale as the recovery of persons with active tuberculosis. The disease progresses, moving from the latent stage to the active one.
This occurs on a time scale that is of the same order as the average lifespan of the carrier of the infection (human).
Tuberculosis can be acquired at random (individual level), that is, as a result of accidental contacts, or through members of an
epidemiologically active cluster that includes at least one actively infected person. The choice of these levels of distribution is
not accidental, it is associated with the observed statistical data of the spread of tuberculosis. Since in the original system there is
no clear separation into fast and slow variables, but the process speeds have different orders, it is necessary to bring the system
(2) to a dimensionless form. To do this, we introduce the following dimensionless variables and parameters:
d S E S1
k t , dt , x1 2 , x2 2 , y1 ,
k k
E1 I k
y2 , y3 , , m , B ,
k k k
where is a small positive parameter. With new variables system (2) has a form:
dx1 x1 x2
d B (1 x1 ) (1 m) y1 n x x ,
1 2
dx2 x 2
(1 m) y2 (1 B) x2 n 2 ,
d x1 x2
dy x x
y1 n 1 2 ,
1
(3)
d x1 x2
dy x2
2 m y1 (1 m) y2 n 2 ,
d x1 x2
dy3 x (1 m) y ,
d 2 3
where y1 , y2 and y3 are the fast variables, while x1 and x2 are the slow variables; x1 corresponds to a population of sensitive
individuals not belonging to the cluster; x2 is a population of latently infected individuals not belonging to the cluster; y1
reflects a population of susceptible people in the cluster; y2 is a latently infected population; y3 is the population of infected
individuals belonging to the cluster. The generated system [4, 5] for (3) is:
dx1 x1 x2
d B (1 x1 ) (1 m) y1 n x x ,
1 2
dx2 x22
(1 m) y2 (1 B) x2 n ,
d x1 x2
x1 x2
0 y1 n , (4)
x1 x 2
x22
0 m y1 (1 m) y2 n ,
x1 x2
0 x2 (1 m) y3 .
The last three equations in (4) determine the zeroth order approximation of the slow manifold (the slow surface) of system
(3). From these equations we have the expressions for y1 , y2 and y3 :
3rd International conference “Information Technology and Nanotechnology 2017” 119
Mathematical Modeling / E. Tropkina, E. Shchepakina
x(t )1 x2 (t )
y1 (t ) n ,
x1 (t ) x2 (t )
x2 (t ) Q x (t ) n x2 (t ) (5)
y2 (t ) n 0 1 ,
x1 (t ) x2 (t ) 1 m
x2 (t )
y3 (t ) .
(1 m)
Since the Jacobi matrix of the fast subsystem of system (3):
g1 g1 g1
y1 y2 y3 1 0 0
g 2 g2 g 2
D m (1 m) 0
y1 y2 y3
g 0 0 (1 m)
g3 g3
3
y
1 y2 y3
has the negative eigenvalues, then the slow surface (5) is stable [4, 5], hence, we can replace the original system (3) by the
reduced one. The reduced system is the projection of the original system on the slow surface (5) with preservation of the
essential qualitative features of the dynamics of the complete system (see, for example, [6-17]).
The reduced system has the form
dx1 x1 x2
d B (1 x1 ) Q0 x x ,
1 2
(6)
dx x x
2
Q0 1 2
(1 B ) x2 ,
d x1 x2
where Q0 n m n ( ) is a number of secondary infections produced by one infected individual in a population each
member of which is susceptible.
System (6) is a homogeneous mixed model in which the infection spread parameter is a function of the parameters:
Q0 n ( ). Recall that the basic reproduction number is defined as c0 Q0 k ( k ). It is easy to see that c0 is the
threshold parameter for the dynamic model (6).
4. Analysis
System (6) has two equilibriums P1 and P2 :
B B B 2 BQ0
P1 (1, 0), P2 , .
Q0 1 (1 B)(Q0 1)
Our goal is to find such conditions on the parameters of the system, reflecting the methods of treatment and preventive
measures, in which the infected and latent infected populations are stabilized at the minimum value. From a mathematical point
of view, this means that the equilibrium of the system is globally asymptotically stable, and their coordinates x1 , x2 should be
as small as possible.
It should be noted that the coordinates of the point P1 correspond to the following situation in real life: all persons are
susceptible to tuberculosis, but not infected (either in active or latent form). In other words, this is the most favorable situation
from the point of view of epidemiology. Therefore, if this point is globally asymptotically stable, then this will be the most
favorable outcome.
The Jacobi matrix of the system (6)
x22 x12
B Q0 Q 0
( x1 x2 ) 2
( x1 x2 ) 2
J
(7)
x22 x12
Q0 (1 B ) Q0
( x1 x2 ) 2 ( x1 x2 ) 2
at the point P1 (1,0) is
B Q0
J P1 .
0 (1 B ) Q0
with the eigenvalues 1 В, 2 1 В Q0 . Thus, according to the reduction principle, the condition Q0 1 В is the
criteria for the globally asymptotic stability of the point P1 . This condition can be written in terms of the main reproductive
number as c0 1.
The Jacobi matrix (7) at the equilibrium point P2 has the form
3rd International conference “Information Technology and Nanotechnology 2017” 120
Mathematical Modeling / E. Tropkina, E. Shchepakina
В 2 В (Q0 2) (Q0 1) 2
(1 В)2
Q0 Q0
J P2 .
(1 В Q0 ) 2 (1 B)(1 В Q0 )
Q0 Q0
For the asymptotic stability of the equilibrium P2 , it is necessary and sufficient that the trace of this matrix be negative and
the determinant positive, i.e.,
В 2 В(Q0 2) (Q0 1) 2 (1 B)(Q0 1 В) 0,
2
В В(Q0 2) (Q0 1) (1 B)(Q0 1 В) (1 В) (Q0 1 В) 0.
2 2 2
The last system can be written as
Q0 1 Q0 0,
(1 В)(Q0 1 В) Q0 1 Q0 0.
Taking into account the physical meaning of the parameters, this condition is equivalent to the inequality Q0 1 B. Thus, we
have the following statement.
Theorem. If c0 1 a disease-free equilibrium P1 (1, 0) (i.e., the point determining the absence of infection in the population)
is globally asymptotically stable. If c0 1, that point P1 (1, 0) is unstable and the equilibrium
B B B 2 BQ0
P2 ,
Q0 1 (1 B )(Q0 1)
is globally asymptotically stable.
It should be noted that in [1] this statement was obtained on the basis of the Hoppenshtadt theorem.
Although the singular point P2 corresponds to the situation when the infected individuals in the population are present, but
for typical values of the parameters, the ratio of the quantity of latently infected individuals to the amount of susceptible
individuals is insignificant. Hence, the situation when this point is asymptotically stable is not so bad from the point of view of
the real situation. In other words, the latently infected individuals will exist but their quantity will not be so high, that is, the
epidemic threshold will not be reached.
5. Correctness of model reduction
The solutions of the original and reduced systems were plotted with the help of Wolfram Mathematica 10.3 and NetBeans
8.0.1. Figures 1-4 show the graphs for systems (3) and (6) with the following parameters values:
1 / 60, 2, 10 5 , n 20, 1.
Fig. 1. The solutions of system (3) with ε=0.00053.
Using a program written in the Java language in the NetBeans 8.0.1 environment, the errors in the deviation of the plots of the
solutions to the complete and reduced systems were determined. For the three cases considered, the errors between the graphs of
the solutions of the original and reduced systems are 0.00010, 0.00012 and 0.00015, respectively. In Figure 5 one can see the
solutions of the original and reduced systems. The figure clearly demonstrates the almost complete coincidence of the solutions
plots, which confirms the correctness of the reduction performed. Thus, the conclusions drawn about the qualitative behavior of
the solutions to the reduced system can be transferred to the original model (3).
3rd International conference “Information Technology and Nanotechnology 2017” 121
Mathematical Modeling / E. Tropkina, E. Shchepakina
Fig. 2. The solutions of system (3) with ε=0.00204.
Fig. 3. The solutions of system (3) with ε=0.00476.
Fig. 4. The solutions of system (6).
6. Conclusion
In the present work, the tuberculosis model has been investigated via methods of qualitative analysis. It has been shown that
the model has several time scales, so it can be represented as a singularly perturbed system of ODEs. The reduction of the
3rd International conference “Information Technology and Nanotechnology 2017” 122
Mathematical Modeling / E. Tropkina, E. Shchepakina
system has been carried out, as a result of which, the original system of five differential equations was replaced by its projection
on the slow integral manifold. It should be noted that, due to the stability of the slow integral manifold, this reduction is correct,
and the reduced system of two differential equations preserves the essential properties of the original model. The conditions
under which the system has the globally asymptotically stable equilibrium have been found. This result means that under the
appropriate selection of treatment and preventive measures, the spread of the tuberculosis epidemic can be completely
suppressed.
Fig. 5. The solutions x1 (t ), x2 (t ) of systems (3) and (6) with ε=0.00053.
Acknowledgements
This study was supported by the Russian Foundation for Basic Research and Samara region (grant 16-41-630529-p) and the
Ministry of Education and Science of the Russian Federation as part of a program of increasing the competitiveness of SSAU in
the period 2013–2020.
References
[1] Song B, Castillo-Chavez C, Aparicio JP. Tuberculosis models with fast and slow dynamics: the role of close and casual contacts. Mathematical Biosciences
2002; 180: 187–205.
[2] Castillo-Chavez C, Song B. Dynamical models of tuberculosis and their applications. Мathematical biosciences and engineering 2004; 361–404.
[3] Vynnycky E, Fine PE. The long-term dynamics of tuberculosis and other diseases with long serial intervals: implications of and for changing reproduction
numbers . Epidemiol Infect 1998; 121(2): 309–324.
[4] Sobolev VA, Shchepakina EA. Model Reduction and Critical Phenomena in Macrokinetics. Moscow: Energoatomizdat Publisher, 2010; 320 p. (in Russian)
[5] Shchepakina E, Sobolev V, Mortell MP. Singular Perturbations: Introduction to System Order Reduction Methods with Applications. Springer Lecture Notes
in Mathematics 2014; 2114: 212 + XIII p.
[6] Strygin VV, Sobolev VA. Effect of geometric and kinetic parameters and energy dissipation on orientation stability of dual-spin satellites. Cosmic Research
1976; 14(3): 331–335.
[7] Gavin C, Pokrovskii A, Prentice M, Sobolev V. Dynamics of a Lotka-Volterra type model with applications to marine phage population dynamics. J. Phys.:
Conf. Series 2006; 55(1): 80–93.
[8] Pokrovskii A, Shchepakina E, Sobolev V. Canard doublet in a Lotka-Volterra type model. J. Phys.: Conf. Series 2008; 138: 012019.
[9] Sazhin SS, Shchepakina E, Sobolev V. Order reduction of a non-Lipschitzian model of monodisperse spray ignition. Mathematical and Computer Modelling
2010; 52(3-4): 529–537.
[10] Pokrovskii A, Rachinskii D, Sobolev V, Zhezherun A. Topological degree in analysis of canard-type trajectories in 3-D systems. Applicable Analysis 2011;
90(7): 1123–1139.
[11] Sobolev VA, Tropkina EA. Asymptotic expansions of slow invariant manifolds and reduction of chemical kinetics models. Computational Mathematics and
Mathematical Physics 2012; 52(1): 75–89.
[12] Korobeinikov A, Archibasov A, Sobolev V. Order reduction for an RNA virus evolution model. Journal Mathematical Biosciences and Engineering 2015;
12(5): 1007–1016.
[13] Archibasov AA, Sobolev VA, Korobeinikov A. Asymptotic expantions of solutions in a singularly perturbed model of virus evolution. Computational
Mathematics and Mathematical Physics 2015; 55(2): 240–250.
[14] Korobeinikov A, Archibasov A, Sobolev V. Multi-scale problem in the model of RNA virus evolution. J. Phys.: Conf. Series 2016; 727(1): 012007.
[15] Sobolev VA. Canards and the effect of apparent disappearance. CEUR Workshop Proceedings 2015; 1490: 190–197.
[16] Lapshova MA, Shchepakina EA. Study of the dynamical model of HIV. CEUR Workshop Proceedings 2016; 1638: 600–609.
[17] Korobeinikov A, Shchepakina E, Sobolev V. Paradox of enrichment and system order reduction: bacteriophages dynamics as case study. Math Med Biol.
2016; 33(3): 359–369.
3rd International conference “Information Technology and Nanotechnology 2017” 123