Sensitivity analysis of the model of tuberculosis

In this contribution, a sufficiently simple, yet realistic model of TB propagation and control is studied in detail. We compute the controlled basic reproduction number R0(uB, uD) for the extended model including the action of two types of controls and analyze the sensitivity of R0(uB, uD) with respect to both controls following the recently devised approach. It is shown that the computed sensitivity coefficients admit a clear epidemiological interpretation and can be used in assessing the efficacy of the respective types of intervention.


Introduction
About one quarter of the world's population is infected with tuberculosis (TB) bacteria. Only a small proportion of those infected will become sick with tuberculosis. However, nowadays TB is one of the top 10 causes of death and the leading cause from a single infectious agent worldwide [1]. People fell ill with TB in all countries and all age groups, but individuals with weakened immune system have a much greater risk of getting infected with TB.
Since decades, there have been long-standing efforts to produce an adequate model of population level dynamics of tuberculosis. We mention the seminal papers [2] and [3] as well as an extensive overview [4]. In our work, we deal with the simple model of tuberculosis introduced in [5]. This model has two control parameters describing the rates of detecting infected individuals from different categories.
Recently, a number of results devoted to the structural analysis of compartmental models were obtained, see, e.g., [6,7]. In particular, the notion of the sensitivity coefficients R 1 was introduced in [6]. The main goal of this study is to apply the results developed in the mentioned paper to the population level dynamical model of tuberculosis and derive conclusions aimed at improving the treatment strategy.
This paper consists of three sections. In Section 2, the structure of the model of tuberculosis epidemiology is considered. In Section 3 we derive the controlled basic reproductive number and perform the sensitivity analysis. Finally, in Section 4 we present numerical results based on real statistical data. differ in some respect. The population is split in a number of subgroups (compartments) that are considered to be homogeneous. That is all individuals within a compartment are identical in their evolution.
When infected with TB bacteria, an individual does not necessarily become sick. Thus, we distinguish two TB related conditions: latent TB infection and (active) TB disease. In the former case the immune system keeps the bacteria under control and so, people with latent TB infection are not ill and cannot transmit the disease. People infected with TB bacteria have a 5-15% lifetime risk of developing active form of TB. If, however, the immune system cannot keep the bacteria under control, the bacteria multiply and destroy tissue. The bacteria usually attack the lungs, but can attack any part of the body. At this point, the person has active TB disease.
We further distinguish two types of active TB: open or smear-positive (SP) and close or smear-negative (SN). In the case of open tuberculosis an inflammation develops inside the lungs. This enables TB bacteria to spread by airborne transmission. In the case of closed tuberculosis, there is no risk of infection as the inflammation is on the lungs periphery and not inside.
Taking into account the above classification, the population is split in 6 groups and the respective states correspond to the number of individuals within the groups: The dynamics of the system include transitions between subgroups as well as the in-and outflows associated with the respective subgroups. There are two control parameters: u D and u B . They describe the rates at which active smear-negative and smear-positive individuals are identified and enrolled in treatment.
The dynamics of disease propagation is described by the following differential equations [5], which are rewritten in order to distinguish between the uncontrolled and the controlled parts: (1) with initial conditions (S (t 0 ), L(t 0 ), D(t 0 ), B(t 0 ), D 0 (t 0 ), B 0 (t 0 )) ∈ R 6 + . Note that we consider the case where there is no influx of individuals from outside into L compartment.
In Table 1, all used parameters are listed and explained. We also show the values of the parameters that were used when performing numerical analysis. These numerical values were obtained from fitting statistical data from Orel Region, Russia to the model (1) and reported in [8].

Sensitivity analysis
We consider the approach to the system analysis based on the concept of basic reproduction number. The basic reproduction number (R 0 ) is defined as the number of infections generated by an infected person in a fully susceptible population. It is a key determinant of the dynamics of an infectious disease. R 0 can be interpreted as threshold parameter: if R 0 < 1 then the disease can not sustain itself in the population and will eventually die out (see, e.g., [9]). That is, to push a disease towards extinction one has to implement a population-level intervention that reduces R 0 as much as possible.
To determine R 0 we present the system (1) in the forṁ where Φ I (x) is the rate at which new infections occur in the population, Φ C (x) represents all flows not related to new infections, and Φ u (x) is the control matrix. The control matrix de-scribes the controlled flows from the infected compartments to other infected compartments (as treated individuals are considered to be infectious, albeit to a lesser extent). Next, we define a disease-free equilibrium. Let x * be an equilibrium state, i.e. Φ(x * ) = 0. Then x * is said to be a disease-free equilibrium (DFE), denoted x DFE , if x * corresponds to the point at which no disease exists. Otherwise, x * is referred to as an endemic infection equilibrium.
We will study the system's dynamics when constant controls (u * D , u * B ) ∈ R 2 + are applied. The linearized model (1) at DFE x = x DFE has the following form: where The matrix A I 11 has non-negative elements. A C 11 is a matrix with non-positive diagonal elements and non-negative off-diagonal elements. The matrices B 11,k have non-positive diagonal and non-negative off-diagonal elements. Furthermore, B 11,k are weakly column diagonally dominant. We introduce the following notation: We are interested in stability of the linearized system (3). The structure matrix of (3) is a block-lower triangular matrix, whose eigenvalues coincide with the eigenvalues of the diagonal blocks. The lower right block is Hurwitz. Thus, the stability of (3)  This spectral radius, parametrized by u D and u B , is called the controlled basic reproduction number: Obviously, we have R 0 (0, 0) = R 0 , where R 0 is basic reproduction number of the uncontrollable system. This approach for finding R 0 (u B , u D ) is based on next generation matrix method [10]. However, in many cases it is difficult to compute an explicit expression for the controlled basic reproduction number R 0 (u B , u D ). Furthermore, in applications it is desirable to be able to compare the contribution of each of the used controls. Knowing that, a policy maker can redistribute the investments in order to focus on a more efficient strategy.
Following the approach developed in [6], we represent R 0 (u * B , u * D ): where R B 1 and R D 1 are the first order sensitivity coefficients, defined as follows: Here, R 0 is the spectral radius of F V −1 and x 0 , y 0 are the right and left eigenvectors of F V −1 corresponding to R 0 . We say that a control u k is efficient if the respective sensitivity coefficient R k 1 < 0. Furthermore, the most efficient strategy is the one with the smallest value of parameter R k 1 , because it contributes at most to the decrease in R 0 (u B , u D ).
Using the formulas above, we compute R 0 , R D 1 and R B 1 for our system (1). The uncontrolled R 0 (0, 0): Note that P is equal to the determinant of the submatrix of V that corresponds to the uncontrollable subsystem or, which is the same, the product of the eigenvalues of the respective submatrix. This value constitutes an important structural characteristics of the system. The expressions for R D 1 and R B 1 are found according to (5). After some algebraic manipulations we write R D 1 as follows: where are the total outflow rates from the compartments L, D, B, D 0 , and B 0 . Note that the structural parameter P 0 can be interpreted in the same way as P. Namely, it is the product of the eigenvalues of the submatrix corresponding to the states L, D 0 , and B 0 . The expression (7) can be analyzed to assess the efficacy of the control u D . Namely, we have that this control is efficient if it holds that P 0 P > β B0 β B k. Similarly, the expression for R B 1 is computed to be This expression has a component similar to that in (7), but also two additional ones proportional to the differences of ratios of parameters related to the treated and untreated compartments.

Numerical computation
To illustrate the theoretical results we first calculated the value of the uncontrolled basic reproduction number, i.e., with u D = u B = 0 according to (4). In doing so we used a particular set of parameter values as in [8] (see Table 1 for the numerical values). The calculations yielded R 0 (0, 0) = 5.1776, which is quite a large value comparable with those of Polio, Mumps and Smallpox. This stresses the importance of choosing appropriate measures to tackle the disease. At the next step we calculated the controlled basic reproduction number R 0 (u B , u D ) for the values of controls taken from the interval [0, 3.5]. The resulting plot is presented in Fig. 1. The value of the basic reproduction number corresponding to the estimated values of controls u * D = 0.2 and u * B = 3.1 is equal to R 0 (u * B , u * D ) = 1.4279 > 1. This value is indicated by the red dot on the plot.
In Figure 1 we see that R 0 (u D , u B ) decreases as the controls grow. The red line separate the region when R 0 is larger then 1 and the region where it is less than 1.
Using the expressions (7) and (8) we computed numerical values for sensitivity coefficients: R D 1 = −8.0227 and R B 1 = −5.9232. We can observe that |R D 1 | > |R B 1 |, which implies that one can achieve the required reduction in R 0 by screening more sputum-negative individuals. Obviously, screening SN individuals requires larger investment. This leads to the formulation of optimal investment problem that is studied in an ongoing project and will be reported elsewhere.